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Combustion—An Aeronautical Science 


HOWARD W. EMMONS* 
Massachusetts Institute of Technology 


SUMMARY 


The basic sciences which are essential to an understanding of 
combustion are briefly presented and their present limitations 
are noted. These limitations imply limitations in our knowledge 
of combustion processes and, thus, indicate areas needing further 
research 


(1) INTRODUCTION 


TT” YEAR 1957 wILL, for all time, be remembered 
as the year in which man took his first steps into 
space. This event, together with the military necessity 
to develop space technology with all speed, gives rise 
to confident predictions of a man in space in the im- 
mediate future, a trip around the moon soon, and 
travel in interplanetary space before the close of this 
century. If all these events come to pass, the twen- 
tieth century will be the one in which man broke from 
his bondage to the earth’s surface, during the first half 
of which he conquered the atmosphere, and during the 
second half he conquered the skies. 

While it is yet too early to predict with more than 
hopeful and wishful thinking the place of atomic energy 
in propulsion, we can say with certainty that the air- 
plane and satellite to date are possible only because of 
the development of propulsion systems built around 
chemical combustion We can also safely 
predict that as long as chemical fuels are our primary 
combustion processes will remain in 


processes. 


source of energy, 
this key position. 

Simultaneously 
aeronautics—and, 
was the growth of engineering science. 
problems posed by aircraft design required such a 
delicate balance of natural forces and material prop- 
erties that only by the merging of science, mathematics, 
and the engineering arts into aeronautical science 
Such fields as dy- 


with the appearance and growth of 
in fact, not trivially related to it 
Many practical 


could working solutions be found. 
namics, aerodynamics, and elasticity received an im- 
mense impetus from the interaction. Combustion 
problems have had the same urgency, have received 
the same careful study from all angles, and have re- 
cently, in fact, been under even greater pressure for 
solution because of their key role in propulsion. In 
spite of this pressure, however, progress has been less 
impressive, less complete. 

I wish to examine the nature of the problems posed 
by combustion processes and to analyze their present 


Presented as the IAS 1958 Minta Martin Aeronautical Lecture, 
M.1.T., March 10, 1988. 

* Jerome Clarke Hunsaker, Professor of Aeronautical Engi- 
neering (Visiting), 1957-1958. McKay, Professor of 
Mechanical Engineering, Harvard University, Cambridge, Mass. 
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status in terms of the factors which currently limit 
man’s ability to understand and predict nature. Thus 
I hope to disclose those combustion areas in which 
future progress is assured and those in which progress 
could be confidently predicted only after certain pres. 
ent day limitations are removed or circumvented. 

COMBUSTION PERFORMANCE 


(2) PROPULSION SYSTEM 


The first powered flight by the Wright brothers in 
December, 1903, owed its success in no small part to 
their development of a powerful, light engine. The 
30 hp. in 210 Ibs., 0.14 hp. per Ib., was an impressive 
development as compared to the other engines of 1903, 
even though it does not sound impressive compared 
to the 3,500 hp. in 3,700 Ibs. or 1 hp. per Ib. in a modern 
piston engine. Similarly, their attained thrust of 
about 90 Ibs.—i.e., 0.43 lb. per Ib. of engine weight 
is no match for the modern jet engine which gives 
10,000 Ibs. thrust, better than 4 Ibs. thrust per Ib. of 
engine weight; or a rocket engine which supplies 
200,000 Ibs. thrust for 4,000 Ibs. of engine, a ratio of 
50 Ibs. thrust per Ib. of engine weight. 

The Wright brothers’ engine, while making a signifi- 
cant development in the use of a combustion process, 
required of the designer little knowledge of the details 
of that process. A gasoline air mixture, if disturbed 
by a spark, exploded, and this was knowledge enough. 
The whole development of the piston engine, both 
spark ignition and diesel, required and actually used 
but little detailed combustion knowledge. The piston 
engine is a masterpiece of machine design through 
strength and vibrations analysis and control, as well as 


materials development and manufacturing process 
control. Only in the knock problem was there a head- 


on encounter with a combustion problem, and, despite 
gallant assaults on all fronts, nature and not man must 
be regarded as the victor. 

In a spark ignition engine, the fuel and air mixture 
as produced by the carburetor automatic flow- 
is compressed in the cylinder and 


an 
metering device 
ignited. The combustion process thus begun propa- 
gates as a flame through the mixture (Fig. 1). If 
the initial temperature and pressure are too high, some 
residual fraction of the combustible will detonate 

i.e., burn at very high speed in a flame which accom- 
panies a shock wave of high amplitude. The transi- 
tion from deflagration to detonation, as this process 1s 
called, remains in the list of unsolved problems, not 
only for the piston engine but for all uses of premixed 
fuels and oxidizers. Not only did nature win the battle 
in this knock process, but also in the much simpler 


process in which nature was kind—namely, the pro- 
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Spark ignition engine combustion (Sloan Automotive 
Laboratory, M.1.T 


portionality between flame speed and piston engine 
speed. No more than qualitatively satisfactory ex- 
planations are available even today. 

I should not imply that nothing useful was learned 
about combustion; quite the contrary. Much em- 
pirical data on composition, temperature, and pressure 
limits of combustion were accumulated. Information 
on fuel additives, fuel chemistry, fuel manufacture, 
Much 


of the thermodynamics of combustion processes and 


and handling has provided essential knowledge. 


products were developed and first used during this 


period. 
With the advent of the jet engine, the combustion 
problem entered a new phase (Fig. 2). At first, the 


combustion problem was assumed to be of secondary 
After all, if one mixed fuel and air, it 
The need to burn more fuel 


importance. 
would burn, would it not? 
at higher speed in less volume soon showed the major 
OS- 


and destructive 


of 


problems. Blowouts, burnouts, 


cillations disclosed large regions our ignorance. 
The growth in performance from heat release rates of 
8 B.t.u. sec.ft.* in steam power plants and 48 B.t.u. 
sec.ft.* in some marine boilers to 12,000 B.t.u. ‘sec.ft.* 
in a modern jet-engine combustor has been accom- 
plished only after many years of development using 
every trick in the book. For reasons that will be made 
clear later, we have not reached that satisfactory state 
in which, by known methods, empirical and analytical, 
and with the assurance of success, we can proceed to 
design a new combustor with specified performance 
characteristics. If the new design differs significantly 
from some tried and true combustor, all the old prob- 
lems of blowout, burnout, and oscillation may again 
return to plague us. 

If the jet engine combustor has resisted satisfactory 
combustion 


too, has the rocket 


These devices have many points 


understanding, so, 
chamber (Fig. 3). 


in common in spite of their differences. Both use a 
liquid fuel. While for the jet engine this is almost al 
ways a hydrocarbon, the rocket immediately used a 
While the jet 


engine combustor must produce a uniform tempera- 


wide variety of exotic compounds. 


ture gas so as to permit the satisfactory function- 
ing of a turbine, the rocket combustor must do so 
After all, the spe- 
per 


in the interest of fuel economy. 


cific impulse defined as the pounds of thrust 
pound per second of ejected propellant carried in the 
vehicle is 2,000 to 3,000 sec. for a turbojet. On the 
other hand, the rocket can currently attain only 200 to 
300 sec. Even future rockets can do little better with 
chemical fuels, unless‘even more exotic species become 
hydrogen-fluorine fuels—or are dis- 


practical—e.g., 


covered—e.g., some magical stabilization of atomic 
hydrogen. Newer fuels can only provide the increased 
specific impulse, and, hence, a smaller and lighter or 
longer range rocket, if the attendant problems of oper- 
ating the combustion chamber at still higher temper 
atures and heat rates can be solved. 

So far, in this brief look at combustion as a propul- 
sion problem, no mention has been made of the solid 
propellant (Fig. 4). From an engine point of view, 
this is simplicity itself. 
enclosed in its case, both storage tank and combustion 
After ignition the burned gases flow through 


The solid propellant grain is 


chamber. 
the grain and out of the nozzle. 
indeed, very attractive and very useful. 
to predict their future. At present the available spe- 
cific impulse is only slightly less than for liquid rocket 
There is a relatively large heavy combustion 


Such propellants are, 
It is impossible 


engines. 
chamber but no machinery of any kind. If problems 
associated with stability of combustion (absence of 
violent oscillations and pressure peaks which cause 
blowups), freedom from danger of detonation, and 
stability of the solid against aging and weather con 


ditions can be attained, then there is a big future for 
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IGNITOR 
COMBUSTION CHAMBER CASE ™ A 
SECTION A-A 
Fic. 4+. Solid propellant combustion chamber 


solid propellants. 


unless the improvement of the specific impulse of solids 
keeps pace with the hoped for improvements in liquids. 

The development of solid propellants with signifi- 
cantly higher specific impulse seems less likely than 
Yet, it would be rash to 


for corresponding liquids. 


predict what frozen free radicals may result from the 
use of materials now under study with, perhaps, the 


literal use of really low-temperature freezing. 


(3) COMBUSTION PHENOMENA 


Let us now turn to a description of what happens, or, 
at least, what we think happens, in a combustion 
chamber. 
bustion chamber fed with nonhypergolic (fluids which 


To be specific, let us consider a rocket com- 


do not ignite spontaneously on contact) liquids. These 
liquids, pumped at about 500 Ibs. per sec. up to per- 
haps S00 Ib./in.* pressure, enter the back of the com- 
bustion chamber through a multitude of small holes. 
Often the liquid jets impinge in pairs. The liquid 
flattens into a disc which then breaks up at its edge into 
Sometimes the liquid is introduced in non- 
In any case, 


drops. 
impinging jets, with or without swirl. 
a shower of fine droplets is produced (Fig. 5).! 
drops of fuel and oxidizer have various histories. Some 
break up in flight by friction with the gas through 
Some evaporate by heat transfer 


These 


which they move 
by convection and conduction from hot products of 
combustion or by radiation from burning drops. Some 
evaporating drops move into gases rich in the opposite 
component and after ignition—an event which raises 
the local mixture temperature until active reaction 
begins—become enveloped in a diffusion flame. Some 
drops collide with other drops and partially break up, 
when of opposite 
composition Oc- 
casional small gaseous pockets of premixed gases are 


and or partially coalesce forming 
some premixed liquid combustible. 


formed either by gas mixing or by evaporation of al- 
ready mixed liquids. These pockets eventually burn 
ina flame. All these events bring the gases to approxi- 
equilibrium the combustion 
This is not complete combustion since 


composition at 


mate 
temperature. 
there is considerable dissociation at high temperature 


This big future will not materialize 
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which readjusts, at least partially, as the gases flo 
through the nozzle. All of this complex of events occy 
for each mass of propellant in the order of a few milf 
seconds. 

While the above description, in the interest of def 
niteness, was given for a particular type of rock 
ccmbustion chamber using two liquid fuels, a som 
what similar series of events occurs in all combustio; 
systems. Essential events are: the introduction of th 
fuel and oxidizer; the breakup and mixing of fuel an 
oxidizer; the ignition of the mixture; the burning , 
the mixture; and the flow of the products of combus 
the chamber. These events may be an 


tion from 


usually are complexly overlapping. Even for solj 
propellants, the events are somewhat the same al 
though the mixing part of the process is initially mor 
or less complete—sometimes the fuel and oxidizer ar 
parts of the same ccmplex molecule—and the breaky 
of the solid, if it does occur, is usually disastrous. 
From an engineering point of view, it may appear 
desirable to try to treat in an overall manner such ; 
This, it 


fact, is what has been rather successfully done o 


complex morass of fine scale phenomena. 


applications of combustion processes throughout the 
Most 
and the bonfire are designed and controlled largely on 


ages. industrial furnaces, household furnaces 


the basis of experience. Provide enough fuel and air 
(an excess of the latter, so as to more completely burn 
the former) for the heating job to be done, provide 
enough space in which the combustion may be com- 
pleted, and provide a sufficiently hot ignition source 
and you are in business. As has happened in so many 
other fields, the aeronautical application puts pressure 
on the designer. In a large rocket combustion cham- 
ber, he must liberate heat by combustion three times 
as fast as a large steam-power station boiler and do it 
in a chamber smaller by a factor 3,000. New difl- 
culties arise and old difficulties take on a new form 
of violence. Thus, it is that engineering which has 
been satisfactorily performed for years with some 
simple overall calculations combined with empirical 
coefficients and safety factors, finds itself unable to 





. 
a 
ow 


Fic. 5. Microflash photographs of formation of spray by two 
impinging jets of water. Views perpendicular and parallel to 
plane of jets are shown (Heidmann and Humphrey ) 
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COMBUSTION 


cope with the new level of precision recuired. Thus 


is born a new aeronautical science. The piston engine 
did not seriously force the issues involved in com- 
bustion. The jet engine did, and the rocket engine 
does even more So. 

Thus. the attainment of 
performance has forced us to look more carefully into 


The vari- 


satisfactory engineering 


the fine structure of combustion processes. 
ous phenomena mentioned above are perhaps best 
organized as in Fig. 6. 

Before discussing the implications and status of the 
overall combustion process, it will be well to recall the 
nature of various parts of the picture. Just as we can 
discuss the combustion chamber as an isolated piece of 
equipment rather than as the entire engine, so we must 
examine combustion phenomena in smaller pieces be- 
fore we can hope to understand the whole combustion 
chamber. This time, however, we will look at separate 
phenomena rather than at separate parts of the device. 


(4) CHEMISTRY 


Clearly, the chemistry is important since without it 
there is no combustion. However, it is amazing how 
little chemistry is needed for many—but not all 
purposes. I do not need to dwell here on the signifi- 
cance and importance of stoichiometry and basic data 
on heats of reaction. Even this part, however, can 
for many purposes—as it, indeed, was for many cen- 
be replaced by simple, ill-stated, and ill-under- 
fire rather than the 


turies 
stood empiricism. The use of 
mortal fear of fire is one of the distinguishing differences 
between men and beasts. Yet, for the major part of 
man's era on earth, he believed either that fire stood 
with earth, air, and water as a major element of our 
universe 01 had no expressed idea at all. Elaboration 
of the fire “‘element’’ idea in the phlogiston theory 
held sway until nearly the year 1S00 when the discovery 
of oxygen and its connection with combustion finally 
made the earlier ideas untenable. 

For overall design purposes, then, only the chemical 
composition of the fuel and oxidizer, the combining 
proportions, and the heats of reaction are needed. 
Even the gross features and some of the smaller scale 
features of the gas movement in combustors can be 
understood by fluid mechanics principles alone, or 
simply by appropriate addition of a heat release at a 
flame front. 

The first additional aspect of chemistry needed in 
the understanding of combustion processes is that of 
the equilibrium composition. Thus, at the high tem- 


perature of combustion products, the equilibrium 
chemical mixture contains not only CO, but CO, Oz, O, 
and, generally, many other compounds containing 
other elements that may be present. In fact, at the 
higher temperatures that would result from more exotic 
fuels, dissociations will be more serious and even ioni- 
zation may play a significant part. The plasma jet 
produced by a high-current electrical discharge is in 


this sense a combustion process because of the im- 
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portant reactions that take place even though a chemi- 
cal reaction is not the primary source of energy. 

The science of thermodynamics developed by Willard 
Gibbs supplies, in principle, all the basic relations 
needed to define and compute the equilibrium com- 
position of any chemical system. For a given temper- 
ature and pressure, the free enthalpy must be a mini- 


mum. Symbolically, 


dG = d(U + pV — TS) < 0 (1) 


for all conceivable processes. For ideal gases—which 
combustion products can quite accurately be assumed 
to be—this equation takes the directly useful form 
of equilibrium constants, expressed in terms of partial 


pressures. For illustration, in the dissociation reaction 


2H,O — 2H: + O, (2) 
the relation between partial pressures at equilibrium is 


Py.’Po, /Puo0? = A(T (3) 
where A(7°) is a function of temperature only and, in 
particular, is independent of the atomic composition 
of the mixture. The equilibrium constant itself is re- 
lated to the thermal properties of the chemical species 
in the reaction through the equation 


din K dT — (f1/ RT?) (4 


where //] is the constant pressure heat of reaction. 
While in principle quantum mechanics describes all 
molecular processes, the calculation of heats of reac- 
tion and equilibrium constants from first (atomic) 
principles is generally too difficult and will continue 
to be so for many years to come for practically im- 
portant reactions. However, empirical data are avail- 
able for reactions encountered in mixtures of the prin- 
cipal elements of most propellant systems— hydrogen, 
oxygen, carbon, and nitrogen. The reason for men- 
tioning quantum mechanics is that some useful results 
are obtained from theory for 10nization and other 
higher temperature processes now becoming important. 
essentially completely 


Equilibrium chemistry is 
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understood for the usual products of combustion. 
There remains in the way of extensive practical use the 
necessity of solving for every temperature and pressure 
of interest a set of equations in the composition con- 
sisting of a linear conservation equation for each ele- 
ment and a nonlinear equation, of which Eq. (3) is a 
sample, for each compound that can be produced 
from these elements. While only simple compounds 
can exist at the high temperature of product gases, still 
there may be a dozen. ‘Thus, practical cases actually 
computed are few in number, but this will probably 
change as computing facilities become more universally 
available. 

The problems of chemical kinetics the 
rate at which reactions proceed—have obvious bearing 
If gases have a residence time in the 


namely, 


on combustion. 
combustion chamber of only a few milliseconds, reac- 
tions which occur in microseconds are instantaneous 
while those taking a second are frozen. Unfortu- 
nately, our knowledge of chemical kinetics is as incom- 
plete as our knowledge of equilibrium thermodynamics 
is complete. The most important reason for this is 
the fact that the overall reaction, as Eq. (2), for ex- 
ample, tells little or nothing about how the reaction 
actually proceeds. For hydrogen and oxygen, for 
example, such reactions as 


OH + H.— H.0 + H 
H + 0.—- 0H +0 | sie 
O+H.—-OH+H 

H +OH+ M-—H.O + ‘| 


are actually found, although it is not yet certain that 
these are correct even for the detonations for which 
they were suggested. 

Thus, during the nonequilibrium, or partial non- 
equilibrium, of a ccmbustion reaction, many fleeting 
intermediate species are present. These species have 
usually not yet been identified. There is no universally 
usable experimental or theoretical technique available 
by which they can be identified. And yet the phys- 
ical, thermal, and transport properties of the mixture 
depend upon their presence. Even mcre important, 
the rate at which the reactions proceed depends upon 
their properties—their molecular properties—their size 
and shape as seen by all the other species present with 
which they interact. 

For each reaction that actually takes place between 
species—as, for example, those listed in Eq. (5) 
kinetic theory, quantum theory, some thermodynamic 
arguments, and, most importantly, experiment agree 
that the reaction rate is given by 

d[A,]/dt = (v."” — vy’) RM [Ay] (6) 
k 
where [4;]| is the moles per volume of species 7; »,”, v,’ 
are stoichiometric coefficients of species s, the former for 
product, the latter for reactant (usually one of them 
is zero): and k is the reaction rate constant of the 
form 
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k = A(T)e~*/*? (7 


where 4, the frequency factor, varies slowly, like 7\: 
and £, the activation energy, is a constant. 

For a given combustion reaction, then, our chemica| 
data would be complete if we knew 

(1) The exact chemical nature of fuel and oxidizer 

(2) The chemical species that will be present among 
the products and during the transient phases of th, 
reaction. 

(3) The physical properties of all species. p, VT 
data from which the corresponding properties of th. 
mixtures may be computed. 

(4) The 
the mixture viscosity, thermal conductivity, and spe- 


transport properties of all species—ie. 
cies diffusion coefficients. 

(5) The thermal properties of all The 
enthalpy and entropy of the mixture must be known 


species. 


in order to follow the state of the mixture. 

(6) The chemical kinetic properties for all possible 
interspecies reactions. This requires identification of 
the important reactions and the determination of the 
activation energy E and frequency factor 4 for each 

For no practically important combustion reactions 
is such extensive, detailed, quantitative data avail- 
able. In the face of this present day limitation in the 
chemistry, combustion processes must be pursued as 
best they can with less complete information. Thus, 
attempts are made to develop an overall reaction- 
rate equation by choosing one best fitted to some 
empirical data. Thus, for example, one finds in the 
literature such global rate equations as 
d[A c]/dt = Kpo,pou,, 1? e7 */*" (8) 


where K = 30.5 X 10” (it.*/sec. mb. °K.) 


Pow Pcatty — [mb. /ft.*] 
E/R = 21,000 °K. 
as the rate equation for the burning of octane. The 
implied chemical reaction is purely fictitious. 

Even if the complete chemical knowledge was avail- 
able and the chemical kinetic problem was, therefore, 
no longer a limitation, the combustion problem would 
still have to cope with the large number of species and 
reactions present. No doubt, some simplifying features 
would be found such as the production of equilibrium 
concentrations of various intermediate species brought 
about by one or a few bottleneck reactions which actu- 
ally control the overall rate. However, there would 
still be a most difficult calculation required in order to 
make use of the information. 

To summarize the knowledge of combustion chem- 
istry, we might say that fuel and product compositions 
are quite well known but all else is in its infancy. 
Thus, severe limitations are imposed on the possibility 
of understanding combustion processes. 


(5) FLurp MECHANICS 


Little emphasis need be placed on the fluid me- 
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chanics problems of combustion for an aeronautical 
audience. All the usual phenomena are present and 
are attacked by the usual experimental and analytical 
procedures. Both laminar and turbulent combustion 
are possible, with the latter playing the major role in 
practical devices. On flame holders, the boundary 
laver is of importance, and perhaps sometimes the 
transition process exerts an influence. Flow separation, 
as always, plays a dominant role when present. 

The gas-flow processes are thus to be treated like 
any familiar gas flow, provided the chemistry permits 
the calculation of the fluid density and other proper- 
ties. The Reynolds Number is still the variable of 
primary importance with secondary roles being played 
by Mach Number, Prandtl Number, ete. 

Even the heat release, as far as the fluid mechanics 
is concerned, is often of secondary importance. Gas 
mixing processes as they occur, for example, in turbu- 
lent diffusion flames are quantitatively identical to 
turbulent mixing without combustion. This equiva- 
lence is so good that highly useful industrial studies 
of the flow through combustion chambers are carried 
out cold, even carried out with water. Thus, eddy 
regions suitable for flame holding, or others, undesir- 
able for flow uniformity, can be located, improved, or 
removed 

As long as tests can be carried out cold, all the 
measurements and instrumentation familiar in a wind 


tunnel are useful. When the hot chambers are tested, 
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flame at one-quarter atmosphere (Fristrom and Westenberg). 
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Fic. 8. Photomicrographs of a low-velocity fuel jet at dif 
ferent distances from the nozzle. Injection pressure, 100 psi 
orifice diameter, 0.020 in; air density, 1 atmosphere (Lee and 
Spencer 


it is much more difficult to get data. Probes must 
be cooled or used for short periods before melting. 
Optical methods are not easy, and sometimes the 
burning gases are not sufficiently transparent. To 
determine the state of the gas, composition measure- 
This requires sampling and subse- 
Recent work has made use 


ments are needed. 
quent chemical analysis. 
of such advanced techniques as spectral analysis or 
mass spectroscopy. While many analyses of com- 
bustion products have been made as a measure of com- 
bustion efficiency, only relatively few measurements 
of the progress of reactions through simple combus- 
tion regions, flames, are available (Fig. 7).* Usually, 
combustion efficiency is measured by the mass aver- 
age temperature which implies total and static pres- 
sure and temperature measurements. 

The one fluid mechanics problem really foreign to 
aeronautical engineers is the impinging fuel and oxi- 
dizer jets. The numerous stability effects of moving 
liquid jets and sheets with viscosity and surface ten- 
sion have been studied analytically and experimentally 
(Fig. 8).* Droplets have been examined, and an ex- 
tensive literature exists on their oscillations, breakup, 
and resultant size distribution (Fig. 9). The overall 
process from smooth jet to cloud of droplets is too com- 
plex to permit prediction of the number, density, and 
size distribution of drops from a knowledge of the 
individual phenomena, but much is known through 
correlation of empirical data using the relatively small 
number of Rey- 
nolds Number, density and viscosity ratios, jet to am- 
bient, and Weber Number (for surface tension). 

Nearly all of this jet information was obtained with- 
As soon as the fuel and oxidizer fog 


necessary dimensionless variables: 


out combustion. 
is ignited, the experimental data are scarce indeed, and 
theories are almost nonexistent. 

The usual limitations of our understanding of fluid 
mechanics—the transition to turbulence, the quanti- 
tative nature of turbulent flows, and the problems of 
flow separation—are all of significance in combustion 
processes. In addition, there are also limitations im- 
posed by the large number of additional fluid-mechan- 
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ics problems arising from the interaction of liquid jets, 
drops, and the simultaneous occurrence of evaporation 


and burning. 


(6) HEAT AND Mass TRANSFER 


The importance of heat and mass transfer varies 
greatly with the specific combustion application. In 
industrial furnaces the heating of a charge may be the 
purpose of the whole operation. In aeronautical com- 
bustion processes the purpose is production of high- 
temperature combustion products, and the engineer 
regards most heat and mass transfer processes as a 
nuisance. Heat and mass transfer means heat losses, 
overheated metal parts, and incomplete combustion. 

The prediction and control of heat and mass transfer 
in aeronautical applications has been attacked _pri- 
marily through boundary-layer studies, both analytical 
sand experimental. The potential usefulness of sweat 
cooling has prompted much of this work. 
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Breakup of water drops in steady stream of air (W. R. 
Land and J. Edwards). 
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Fic. 10. Composite heat and mass transfer correlations for 


spheres (Ranz). 


Recently, of course, all these studies have been use- 


ful in the re-entry problem. Here we have a ‘‘com- 


bustion”’ problem of really manageable proportions. 
The high-temperature gas behind the bow shock dis- 
sociates and recombines, transfers heat to the re-entry 
The 
chemical reactions and ionizations are so few in number 
that it will be possible to follow them in detail. The 
fluid motion is in large part laminar, thus making 


body within a shock and boundary-layer region. 


analytical treatment hopeful. The relative simplicity 
of this problem and the widespread effort now being 
made will within a few years yield satisfactory solu- 
tions. No such optimism can be justified for the gen- 
eral combustion problem. 

While 


processes are not wanted, they are essential as far 


most aeronautical heat and mass _ transfer 
as the new charge of fuel and oxidizer is concerned. 
Only through these transfer processes can the new 
charge be mixed, evaporated, and ignited. These 
processes are like the fluid mechanics of jets and drops, 
often understandable as small separate pieces, they 
become impossible when the complex is viewed as a 
whole. 

The progress that has been made with complex 
industrial heat and mass transfer processes, largely 
by chemical engineers, has been achieved through the 
accumulation of carefully taken experimental data (1 
will not term it accurate, however, since sometimes 
experimental difficulty makes standards low) corre- 
lated with 
(Fig. 10).° 


the analysis of data has been the introduction of the 


the appropriate dimensionless variables 
A key step in the measurements and in 
idea of a thermal circuit. Heat passing to a wall, 
through the wall, and then away from the wall on the 
other side constitutes a thermal circuit. By inde- 
pendent study of each of the parts of such circuits, 
their thermal resistance characteristics (usually ex- 


pressed as conductance—1.e., heat-transfer coefficients 


can be quantitatively handled. This process has 
made reasonable engineering design possible with 


phenomena for which little or no analytical under- 
standing has yet been achieved. The heat transfer 
within combustion chambers is so complex that little 
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has been done even with the accumulation of imprecise 
data, and almost no correlations have been attempted. 


(7) RADIATION 


The role of radiation in combustion processes is in 
an uncertain state. Radiation is always produced 
during the chemical reactions and in the hot products of 
combustion. In a sense, it 7s one of the products of 
combustion, a shower of photons being produced along 
with the various molecular species. The gas radiation 
arises from two sources—the excited states in which 
molecular or atomic species first are produced, and by 
thermal radiation of the final product gases at their 
high temperature. 

Radiation is used through spectral analysis to study 
the chemical species present. Of greatest importance, 
however, is the radiant heat transfer. Radiation to 
metal parts is a thermal circuit, parallel to the con- 
vection heat transfer. While the energy transferred 
by radiation is not a large fraction of the heating 
value of the fuel (about 5 per cent of the heat of reac- 
tion appears as radiation at the boundaries of the com- 
bustion chamber), it is not insignificant. In one re- 
spect, at least, it is of major importance 
Unless thermocouples or 


namely, in 
temperature measurements. 
other sensing units can be shielded from the radiation, 
a serious error may be made in the indicated temper- 
ature. In fact, in a combustion process even the sign 
of the error is not certain since the flames radiate vari- 
ous spectral bands at very high temperature (some 
well above the flame temperature itself) while cool 
walls radiate at a relatively very low temperature. 

The origin of radiation goes well into the molecular 
and atomic structure of While much is 
known about energy levels in molecular species, this 


matter. 


knowledge is useful in combustion only in chemical 
analysis and not in overall energy and flow consider- 
ations. In fact, the role of radiation in flame stabili- 
zation still presents open questions. 

Leaving aside the nature of the origin of radiation as 
too complex for direct incorporation into performance 
calculations, we turn now to the properties of the 
radiation itself. Electromagnetic radiation in a com- 
bustor has all the usual properties: wavelength and 
wavelength dependence of intensity and polarization, 
and interaction with matter through wavelength de- 
pendent refraction, reflection, transmission, absorption, 
and emission. Such a complex of phenomena again 
makes a complete understanding hopelessly difficult. 
The final blow to any hope for reasonably accurate 
measurements or calculation is the fact that radiation, 
absorption, and reflection by most metals occur in the 
first few atomic layers, and, hence, minor chemical 
or physical changes of the surfaces present make for 
large changes in the results. 

Under such adverse conditions, it is not surprising 
that little has been done with the measurement and 
analysis of radiation in combustors. The progress 
made in combustion radiation studies of furnaces by 


chemical engineers makes it clear that more can be 
done when By assuming grey radiation 


and grey properties (‘‘grey’’ means a fraction of a black 


necessary. 


body radiation independent of frequency) and _ per- 


forming such integrations as 


. . 
a = (m/rl | cos O(¢ 
JAdJI 


the radiant heat transfer from a volume of gas to a 


r*) dV d.l (9 


surface area .! can be calculated as 
QO = amoTgil (10 


where the notation is shown in Fig. 11; -1 is the area 
of interest, m is the gas absorption coefficient, V is the 
entire volume of gas seen from d.1 (therefore a complex 
function of .1). 

While correlation techniques are indirectly implied 
by engineering radiation heat-tiansfer calculations, 
they have played little direct part in radiation studies. 
In a sense, the evaluation of angle factors and other 
integrals like Eq. (9) by experiments on models makes 
use of correlation techniques for black body radiation. 

Almost no attempt to make exact calculations of gas 
flow fields with radiation has been made outside the 
field of astronomy. Even the simplest case of radiative 
equilibrium leads to difficult integral equations. The 
modest magnitude of the radiative energy transfer 
in combustors and the large magnitude of inherent 
errors because of uncertain data have made the ex- 
amination of these difficult integral equations unat- 


tractive. 
(S) COMBUSTION 


We now have examined briefly the various fields of 
knowledge which come together in combustion. In 
no practical combustion problem is it possible to sim- 
plify the phenomena to a single one. Thus, there is 
always present a degree of coniplexity beyond what 
can be readily handled by simple techniques—either 
analytical or experimental. 

The general problem, from an analytical stand- 
Although 


point, is shown in Eqs. (11-19) (see p. 738). 
kinetic 


these equations are frequently derived by 
theory,’ such a procedure is not necessary since they 





are essentially continuum equations.’ By proper 
selection of terms and boundary conditions, these 
Po 
PY 
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equations embrace practically all problems involving 
fluid continuum, except for the recently popular sub- 


ject of magnetic and electric effects. 
Only the simplest combustion problems have been 


examined in so exact a way as to use all the obviously 
Even then, there 


necessary terms of these equations. 


are usually additional simplifying assumptions added. 
To illustrate, let us look at a few simple problems. 

Consider first the calculation of flame speed—the 
velocity at which a reaction front moves through a 


premixed mixture of fuel and oxidizer. 


A presumably 


very simple chemical case is the hydrazine decom- 


position reaction 


NoH, —! 


2NH2 


(20) 


for which the thermal, transport, and kinetic data are 
After considerable simpli- 


presumed to be complete. 


fication we are led to the two coupled equations 


de,/d0 = —[(aA),/s*] X 
} [Xi exp (—@/0)|/[(@ — 1) + (1 — @)al} (21) 
dX\/d@ = Le) (2X1 — « — Xie) + 
[(@—1) + (1 — aly (22) 
where 
€: = flowing mass fraction of hydrazine 
a = thermal diffusivity 
A = frequency factor of Eq. (7) (evaluated at the 
final temperature, indicated by subscript /) 
S = flame speed 
X, = mole fraction of hydrazine 
6 = (7T/T,) = temperature ratio to final temper- 
ature 
6, = (E/RT;,) = activation temperature ratio 
6 = initial mixture temperature ratio 
Le = Lewis Number = a/D 
D = diffusion coefficient of NH»in N.H. 


To satisfy the initial and final boundary conditions 
a specific value, an eigenvalue, of the flame speed, § 
must be chosen, provided the reaction rate is negligibly 
slow at the mixture temperature. Calculations show 


that for various Lewis Numbers we get results as fol- 


lows® for hydrazine originally at 423°K.: 
Le 0.444 0.889 co Experimental 
S (cm./sec.) 27.9 39.3 95 200 


The poor agreement for all possible diffusion coef- 
ficients suggests that the simple chemical reaction 
mechanism assumed is probably incorrect.* 

A few other reactions have been thus studied. Some, 
like the ozone decomposition flame, seem to check quite 
well. For most combustion flames, however, there is 
too little chemical and physical data even to attempt 4 
solution. At 
flame speed is the only practical approach to evalu- 


present, experimental measurement of 


ation of overall combustion mechanisms. 


Next, consider a diffusion flame, as in Fig. 12,° 


where 
fuel is supplied through the plate as it is evaporated 
from the surface. The Eqs. (11-19) written for a 
boundary layer still look rather formidable. How- 
ever, the complete solution can be constructed if it 1s 
assumed that 

(1) The release of heat occurs proportionally to the 
rate of destruction of fuel and oxidizer and the rate 
of production of combustion products. 

(2) Prandtl Number = Schmidt Number = 1. 

(3) pu = constant. 
The resultant rate of combustion and drag on the 
plate are shown in Fig. 12. The combustion coefficient 
C, is defined as the mass rate of burning of fuel divided 
by the mass rate of flow of free-stream air through 


* Recent work introducing a chain reaction mechanism in 


place of Eq. (20) gives excellent agreement with the experimental 


results. !? 13 
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a normal area equal to the plate area, while the heat 
ratio, B. is the enthalpy of a unit mass of oxidizer 
both thermal and chemical) relative to the fuel sur- 
face divided by the energy required to produce a unit 
mass of fuel vapor. Thus, with the above assumptions, 
no reference need be made to the actual reaction rates, 
provided that only combustion rate and drag are de 
sired. This means that for some aspects of a combus- 
tion problem the detailed chemical information is of 
little or no importance. In the diffusion flame, the 
transport properties are all important. However, if 
we want to know the spacial distribution of species or 
temperature, or the thickness of the burning boundary 
layer, reaction rate data are essential. 

Similar calculations to the above are suitable for 
many diffusion flames, such as burning drops. 

As a final example, I want to mention the missile 
re-entry problem. Here, as indicated before, the bow 
shock wave on a body covers the body with flowing, 
reassociating molecular fragments. The shock itself 
already produces the equilibrium distribution of frag- 
ments, atoms, radicals, molecules, ions, ete. Through 
the inclusion of the appropriate terms of Eqs. (11)—(19), 
the solution for the whole flow can be found numerically. 

By an assumed pressure distribution on the body 
start with Newtonian), the compressible Bernoulli 
equation and the thermodynamic, and, if important, 
chemical kinetic laws permit the evaluation of state 
and velocity along the body (as a streamline). Fol 
lowing this, continuity provides the distance normal to 
the body to the “‘next’”’ streamline and, hence, makes 





tal > eo . . : 
Fic. 12. Combustion in a boundary layer by evaporation of 


fuel from a plate 
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Fic. 13. Missile re-entry flow (test calculation 


it possible to construct a second streamline a small dis- 
tance removed from the body. Next, the momentum 


equation normal to the streamlines, 


Op On —pV?/R (2: 


makes possible the calculation of the pressure increment 
between streamlines in terms of the known state, p, Ss, 
velocity, V, and streamline curvature radius, RX, and 
thus the pressure distribution on the new streamline is 
determined. Finally, we start the process over again 
to construct the entire flow streamline by streamline. 
By matching the streamlines thus constructed to the 
free stream, a shock position is found, but, of course, 
the shock pressure may not agree with the pressure cal- 
culated by the construction. Thus, it is that a correc- 
tion to the initial body pressure distribution is found 
which permits the whole calculation to be repeated. 

The success of this process is shown by Fig. 13 where 
an incorrect cone surface pressure distribution was 
assumed, and after two iterations the solution had 
settled down to the correct one—the exact solution 
being known in this simple case. 

The method of solution described here“ 


is one of a 
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Fic. 14. 
correlation with Damkohler first number. 


number that have recently appeared. It has the ad- 
vantage of working without any restrictive assump- 
tions on the properties of the gas and, hence, permits, 
with a minimum of additional computation, the use of 
real gas properties, including chemical reaction rate 
data as they become available. This re-entry problem 
promises to be one of the few wholly solvable com- 
bustion problems. 

Because of the inherent difficulty of a direct analyt- 
ical attack, combustion, as it occurs in combustion 
chambers, will be studied and developed for many 
years to come by the test of trial and error. Thus, 
the use of models and a study of model laws will always 
be in the foreground of test evaluation. 

The differential equations themselves suggest a 
number of essential parameters. 


Reduced frequency wl 

Reynolds Number Re = VLp yu 

Froude Number Fr = V/(Lg)!/? 

Schmidt Number Sc = uw pD| Lewis Number Le = 

Prandtl Number Pr = C,u/kf Se’ Pr 

Isentropic exponent y = C,/C, 

Mach Number JJ = Va 

Molecular weights for species s, .1/, 

Stoichiometric coefficients for species s in reaction 
r, as a reactant »,,’, as a product 7,” 

The frequency factor for reaction r, .1,(7)) 

The activaticn-temperature ratio for reaction 7, 
Or = E/RT 

Damkohler’s second number for reaction 7, Daz, = 
Of/t,f 

Radiation number RK = eo7* C,pu (o is Stephan 
Boltzmann radiation constant) 


Many of the above dimensionless variables are func- 
tions of the state, so that one must further define at 
what point in the combustor the number is to be 
evaluated. 
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Combustion efficiency for turbojet combustion chamber 
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In addition to this list of S + s + 2rs + 3r dimen. 
sionless variables arising from the differential equations 
there are also others which arise from boundary cond. 


tions. 


Weber Number We = V°*L,y (for liquid gas inter 
face) 

Initial composition ]Y,, (fuel air ratio) 

Evaporation number B = Q X 

Various geometric ratios 

Various property ratios for fuel and oxidizer 


The great simplification, usually brought about by 
dimensional analysis and indicating clearly the vari- 
ables that are really important, is also connected with 
the fact that the independent variables of a problem 
are reduced to a manageable number—say, from five 
to two. 

In combustion the general problem is not improved 
much because of the initial very large number of 
variables and their close connection with the nearly 
unknown chemical mechanisms. 

The only way now known to cope with this situation 
is to find legitimate (or, sometimes, other) reasons for 
splitting the problem into smaller pieces, each of which 
is simpler. Sometimes this split can be a physical one 
as, for example, the separate consideration (partially 
legitimate) of the liquid jets and their breakup, the 
burning of fuel droplets, and the flow of reaction prod- 
ucts through the nozzle. Sometimes this split is simply 
an attempt to select the one (or a few) phenomenon 
which exerts the major control on the overall process 
The chemical species and reaction-rate variables are 
“eliminated”’ from consideration by, first, using only 
one fuel and oxidizer—sometimes at a fixed fuel-air ratio 
and initial temperature——and, second, by assuming 
the existence of a ‘‘bottleneck’’ reaction, specified 
adequately as a global reaction. All the chemical 
complexity is thus replaced by a single Damkohler first 
number, defined as the ratio of fluid residence time 
L/V, toa characteristic chemical reaction time, Tohem 


Da, L inom (24 


While this is the definition, it is never used because 
the chemical time Tehe,. is neither well defined nor ever 
known. Rather, it is assumed that the flame speed or 
some other performance characteristic is a measure 
of the chemical times, and, by a study of experimental 
results, the dependence of the chosen characteristic on 
pressure is determined in the form p”. Thus, as actu- 
ally used, the Damkohler first number is expressed 


in the not so dimensionless form 
Dua = Lp"/V (25 


A variety of “‘most significant’’ variables has been 
produced in this same approximate way, and _ these 
have shown in the hands of their authors some merits 
However, to date 


as correlation variables (Fig. 14)."! 
the amount of experimental data on geometrically simi- 
lar burners under a range of operating conditions is 
net at all adequate to really explore the effectiveness 
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COMBUSTION 


and limitations of the proposed correlations. Essen- 
tially, it is still necessary to develop combustion 
equipment in prototype for every application. 

Progress with combustion processes, both in _per- 
formance and understanding, has been great, especially 
in the past 15 or 20 years. As the fields of fluid me- 
chanics, heat transfer, thermodynamics, and reaction 
kinetics have moved ahead, the application to com- 
bustion has followed. Similarly, the limitations of 
these separate fields immediately imply limitations in 
the corresponding part of combustion processes, be it 
experimental or theoretical. Even when the separate 
fields in a given combustion process are well within the 
domain of our understanding, the latter may be out of 
bounds because of the complex interactions between 
fields. 

Thus, the essential limitations of combustion, in a 
sense, are the same as in all engineering applications 
many inde 


of classical physics and chemistry: too 


pendent variables. Independent variables, up to five 
or six, make for an easy problem. <A.n infinity of inde- 
pendent variables of a single or several types makes 
Kinetic theory of 


3ut what can we do 


statistical methods applicable. 
gases is the classical example. 
thousand different and 
To 
The engineer's problems are 


when there are ten to a few 


individually important independent variables? 
date there is no answer. 
almost all of this kind, and his method of sclution 1s 
to use “engineering judgment.» The business manager 
uses “business judgment.’ The politician uses “‘politi 


cal acumen.”’ The military uses “‘generalship,’’ per 
haps supported by operations research. The numerical 
Monte Carlo methods. On_ the 


scientific edge of some of these branches of man’s at- 


analyst may try 
tempts to control nature (including himself and _ his 
fellow men), there is a growing use of the large-scale 
computer. On the empirical edge of these branches 
the process is “‘geod judgment’’; the name given to the 
process of selection of the right course of action, as 
determined a posteriori. 

the aeronautical stands at 


Combustion science 


318, JPL, 
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that point in its history where engineering judgment 
can expect some, but not very extensive, help from 
the systematic, rational, and logical processes known 
as science, and we can look forward to the day, though 
distant, when our knowledge of the handling of com- 
plex systems, and of this field in particular, will be 


sufficient for the rational design of combustion devices. 
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Approximate Three-Dimensional Solutions for 
Transient Temperature Distribution 
in Shells of Revolution 
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SUMMARY 


An approximate three-dimensional heat conduction equation 
is formulated for thin shells of revolution by making assumptions 
similar to those of the theory of thin elastic shells. It is then 
shown that the resulting partial differential equation is separable 
for all shell shapes and that solutions may be obtained in terms 
of well-known functions. It is also possible, for a broad class 


of problems, to solutions as combinations of a one- 


dimensional slab solution and suitably defined correction func- 


express 
tions 


INTRODUCTION 


ppg OF the rapid developments in the area of 
high-speed flight, problems of thermal stress and 
deformation in aircraft structures have, in recent years, 
received a great deal of attention. The calculation of 
these thermal stresses and deformations requires a 
knowledge of the transient temperature distribution 
through the structure which may easily be obtained 
for simple elements such as slabs, rods, and _ plates. 
In addition, useful approximate methods have de- 
veloped, and, in particular, Biot has brought to bear 
on this problem the powerful methods of the calculus 
of variations. ! 

The calculation of transient temperatures in shells 
of revolution is in many cases quite difficult because 
the heat equation is separable only in a limited number 
of coordinate systems.” It therefore seems worth 
while to investigate whether a simpler form of heat 
equation may be obtained by making assumptions 
similar to those of the theory of shells. In this manner, 
we obtain a formulation of the heat conduction problem 
in which the approximations are of the same order as 
those of thin shell theory, which is used subsequently 
in determining thermal stresses and deformations. 


(1) THe Heat CoNnpucTION EQUATION FOR SHELLS 
OF REVOLUTION 


Referring to Fig. 1, and using the cylindrical co- 
6, z, the parametric equations for the middle 


et 


ordinates /, 
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surface of the shell may be written as 


r= = 2(&) | 


r(é) and g 


where & is any convenient parameter. We now con- 
sider the triply orthogonal right-handed coordinate 
system & 6, ¢ where ¢ is measured along the outer 
normal to the middle surface. With respect to this 
system, the square of the element of arc may be written 


ds? = a®{1 + (¢ r,) |*dé? + r?[1 + (¢/r¢) ]*d6? + de 


(? 
where 7; and 7, are the principal radii of curvature of 
the middle surface and 

a = [r’? + 2/2]1/2 (3 


while the primes denote differentiation with respect 
to &. 
The principal radii of curvature are given by 


a/g', rp = r/sin ¢ (4 


where ¢ is the inclination of the normal to the middle 
surface. 
From Eq. (2), we obtain the covariant metric tensor 


associated with the (, 6, ¢) coordinate system as 


a?(1 + ¢/r;)? 0 0 
4 = 0 r\1 + C/r)? O (5 
QO 0 l 


and the contravariant metric is, therefore, 














Fic. 1. 
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ii : ea 
0 0 
el + t/Pr.)* 
= l (6) 
0 0 
r*(1 + ¢/t4)* 

8) 7] l a 

In terms of general curvilinear coordinates x‘(7 = 1, 2, 


}), the heat conduction equation for a thermally 


” 


jsotropic medium with constant thermal properties 


may be written as 


V?7T = (1/a)(OT/Ot) (7) 


where a is the thermal diffusivity. The Laplacian 
operator V°7" is defined as 


V? = g"™T im (8) 
where 7°,,,, denotes the second covariant derivative of 
T with respect to x’ and x”. Since 7 is a scalar, 


its first covariant derivatives are the partial derivatives, 
and we have 


T,, = OT /Ox' (9) 


Tim may now be obtained by covariant differentiation 


| jo°T 
ae(1+E/r_)? LOE 


For the of thin shells, the wall thickness is 
assumed small compared to the smallest of the principal 
and ¢ 


In addi- 


case 


radii of curvature, and the terms ¢/7; rg are 
therefore negligible in comparison to unity. 


tion, we note that 


(¢/ra)¢’ cot ¢ = (£/r,) [In (sin ¢)]’ < (r’/r) 
and, similarly 
(¢/rg)(e"/e’) = (¢/rg) [In 9’)! < (a’/a) 


The heat conduction equation for thin shells is now 
obtained by making use of these characteristics and 
neglecting appropriate terms in the coefficient functions. 
The result is 

1 /a®) {(0°7/0e?) + [(r’/r) — (a’/a)](OT/0E)} + 
(1/r?) (027/002) + (027/d¢2) + [(1/r_) + 


(1/r5)] (OT /0¢) = (1/a)(OT/dt) (15) 


In addition to the above simplifications, it seems reason- 
abie to expect that the term [(1/r;) + (1/74) ](O7/0¢) 
will not have a strong influence on the variation of the 
temperature with ¢. This term represents the effect 
of area variations in the ¢ direction on the heat balance, 
and since the wall thickness is assumed small, we 
would expect the term [(1/7;) + (1/7r¢)](O7/O¢) to 
be small compared to 0°7/0¢°. To investigate this 
further, we consider a case where the temperature 
does not depend on ¢ and 6, and we define new dimen- 


74 


> 


of the covariant vector 7,,. We then have 
re o*T js | oT 
Tim = —i.¢ (10) 
Ox'Oxn” tim § ax 
, ps ; ae . ; 
where Yim ¢ denotes the Christoffel symbol of the 
m : 
second kind. The Christoffel symbols of the first 
kind are defined as 
(77, a] = (1/2) [(Og;n/Ox’) + (Og;,,/O0x') — 
(Og;;/Ox")] (11) 
and the symbols of the second kind are 
st ; 
,..¢ = gl, (12) 
\ ij ( g°"[7j, ] 
For an orthogonal coordinate system with g' 0 for 
i ~ j, Eq. (12) reduces to 
) ii| = g**(7j, s] (no sum on s) (13) 
Eqs. (5), (6), (11), and (13) may now be used to 


calculate the Christoffel symbols for our coordinate 
system. The results can then be used in Eqs. (7), 
(S), and (10) to obtain the heat equation in coordinates 
of revolution. This equation may be written 
7 . o°T 
O6- oc 


l/r; ‘s 1/ro \¢ 107 (14 
+0 /%%, 1+¢/rejO¢ a Ot 
sionless variables as 
t’=ft/kh and ?’ = at/k* 
Eq. (15) then reduces to 
(0?7/O¢"?) + [(h/r~g) + (h/r) (OT /0¢") OT /Ot’) 


Because of the presence of the term [(//r;) + (h/r)], 
it is clear that the second term on the left-hand side 
of this equation will be negligible compared to the 
first, as long as O7'/O¢ does not become much larger 
than 0°7/0¢°. Since such a situation is unlikely for 
a thin shell, we conclude that the term in question can 
be dropped in Eq. (15) and the heat equation for thin 


shells becomes 
a) (OT /0E)} + 
(1/a)(OT/Ot) 


(1 a?) }(0?T O&) + [(r’ r) — (a’ 


(1/r?)(0?7T'/00?) + (0°T/O¢?) = (16) 


The accuracy of this last approximation is investi- 
gated for specific cases in Appendix (A). In expressing 
the boundary conditions for typical problems, we will 
require the physical components of the heat flux 
The covariant components of this vector are 


qi = kT,, = —k(OT/Ox*) 


vector. 


and the physical components are obtained by multiply- 
ing the covariant components by the appropriate 


metric coefficients; thus, we have 
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dt = V oan, de = V 92 qe, qr = V g54q, 


Substituting the appropriate values of q; and g’ in 
the above expressions and introducing the shell ap- 
proximations, we obtain 


ge = —(k/a)(OT/0E) 
do = —(k/r)(O0T/08) 
qc = —k(oT O¢) \ 
APPROXIMATE 


PROPERTIES OF THE 


EQUATION 


(2) GENERAL 


With the simplifications introduced above, the tran- 
sient temperature distribution in a shell of revolution 
is a solution of Eq. (16), together with appropriate 
boundary and initial conditions. Since the tempera- 
ture function must be periodic in 6, we can express the 
solution of Eq. (15) in the form 


x 


> V,(é, ¢ 0) cos n0 + 


n=0 


T(é, ¢, 0, t) = 


zz Vn»(é, ¢, t) sin n@ (18) 
m= 1 
where the two series represent the even and odd parts 
of the temperature function. Since the determination 
of the functions |’, and V,,, involves identical problems, 
we will consider hereafter only the case of even tempera- 
ture which includes the axially symmetric problem 
as a particular case for = 0. We therefore substitute 
in Eq. (16) the first series on the right-hand side of 
Eq. (16), and we obtain 


1/a?{(0?V,/d8) + [(r’/r) — (a’/a) (OV,/d&)} — 


(n?/r?)V_ + (0°V,/08*) = (1/a)(OV,/ot) (19) 


For most cases of interest in aircraft structures, the 
boundary conditions on 7 may be expressed in the form, 
Ciqy + GT = flé, 


= 2( 
Cs4¢ “bh Ci = \ )) 


| 

S 

wr D> 
II 

wr 


It is seen that the boundary conditions (20) allow 
prescribing an arbitrary or convective heat flux or a 
given temperature on the inner and outer surfaces of 
the shell (¢ = +//2). 

Eq. (20) may easily be expressed in terms of J’, if 
the function f is expanded in terms of cos 76 or sin 6; 


thus, we obtain 


—CiR(OV»/O5) + C2Vn = frlé, 0), § = +h/2Y (91) 
¢ (2 


—C;(k a)(OV, O&) + Cilia = (), é = £1, &: 


(1. *)f 7 


From the first of Eqs. (21), we see that, in general, the 
boundary conditions on J’, will be nonhomogeneous. 
Such nonhomogeneous conditions may easily be re- 
moved by taking I’, in the form 


where /,(&, ¢) = f(E, 0, t) cos nédé 


VE, £,8) = Alé, £, 0 + 2e©fnlé, t) (22) 


where the function g(¢) is chosen such that the boundary 
conditions on H, are homogeneous; 


this procedure 


SPACE 


1958 


SCIENCES DECEMBER, 
will be illustrated later. When Eq. (22) is substituted 
into Eq. (19), we obtain a nonhomogeneous equatio; 
for /7,, thus 


(1 a?) }(O7/7, of) + [(r’/r) — (a’ a) |(O7T,,/0&){ — 
(n?/r*) HT, + (O71T,,/067) = (1/a)(OlT,,/0t) + Q,(€, ¢,t 


(99 


-) 


It is useful, at this point, to investigate the existeng 
of separable solutions for the homogeneous part of 
Eq. (23). 
used to obtain the particular solution appropriate t 
We, therefore, consider the equation 


Such solutions may then be conveniently 


any given Q,. 


(1 a?) } (O° O€) + [(r’/r) — (a’/a) JOH O§} = 
(n?/r*)H + (0°7H7/0¢?) = (1/a)(OH/ot) (24 


where the subscript 7 has been dropped for convenience 
By virtue of our choice for the function g(¢), the 
boundary conditions on // will be of the form 


—Ck(0H/0¢) + C.7 = 0, € = +h/2 \ a 
—C;(k/a)(0H/dt) + GH =0, ¢=&,6f 
Assuming solutions of the form /7(é, ¢, t) = ®(&) x 


Z(¢)F(t) and substituting in Eq. (24), we obtain for 
®, Z, and F the following three differential equations: 
a) |b’(é) + 

a?[\2 — (n?/r?)]@(E) = O (26 


e’(é) + [(r’/r) — (a’ 


Z"(¢) + y*Z(¢) = 0 (27 


F’(t) + a(X? + 7?) F(t) = 0 (28 


where \” and y” are eigenvalues to be determined from 
the boundary conditions. The boundary conditions 
on #(é) and Z(¢) are obtained from Eqs. (25); thus, 


we have 


—CkRZ"(6) + O25) = 0, § = th/2 | 
—Ci(k/a)b'(E) + Cyb(~E) = 0, E = &, kof 


(29 


It is now clear that separated solutions of Eq. (24) are 
possible for all shells of revolution. This is a sub- 
stantial advantage over the exact heat conduction 
equation which is separable only in a few coordinate 
systems. Thus, we now have a convenient method of 
obtaining temperature solutions which are approximate 
to the same degree as the equations of the theory of 
shells used in the determination of thermal stresses 
and deformations. 

There will, in general, exist an infinity of values 
A»; and y, for which the boundary conditions are 
satisfied. 
time function F will be 


As a result, we see from Eq. (28), that the 


| exp [—a(A,;7 + 7,7)é] (30 


So that the homogeneous solutions for //7,(é, ¢, t) will 
be of the form 


H,(é,¢,t1) = dS Anspexp [—a(Ans? + yp)t] X 
s=1 1 


p= 
] 


®,(€)27,(6) (51 


ns\S/*™ PX) 


The coefficients A,,, must therefore be determined 
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It is therefore of 


irom the initial conditions on H. 
which 


interest to determine the conditions under 
Eq. (26) will give a set of orthogonal functions as 
solutions. Multiplying Eq. (26) by v/a and rearrang- 
ing terms, we obtain for a given characteristic root A 


32) 


d dé) \(r a)®,’(E)] + [ard.? — n?(a/r)]&(E) = 0 


It is clear that Eq. (82) and the boundary conditions 
299) form a Sturm-Liouville system, and we can 
therefore state* that any two solutions ©, and ®, of 
this system will be orthogonal with respect to the 
weighting function ar as long as the functions ar, 
ra,and @/r are continuous in the interval & < & < &. 
In addition, if a/r does not change sign in the interval, 
every characteristic root \, will be real. It should be 
noted that, since 7/ a must be continuous in the interval, 
ar cannot have roots in this region, and therefore the 
condition for real characteristic numbers is automati- 
cally satisfied. The proof of these properties is well 
known and will not be repeated here; we will, however, 


make use below of the expression 


ar — Ay) | ar® Pdi = 


| Vr a)( P,P ! == ©, '@ )} (33) 
which may easily be obtained from Eq. (26). 
We will need later to determine the norm of the 
functions ®, given by 
N, = ar®-dé (34) 


This can be done by making use of Eq. (33) as follows: 
consider two functions ®, and ®, obtained by choosing 


neighboring values of \ such that A; = A, + 6A,. We 


then have 


®, = &, + (d&,/dd,)bA,, 


and A? — A,? = + ddA,” 


2,6 


Introducing these values in the left-hand side of 


Eq. (33), we obtain for the left-hand side the expression 


2r,6A, + 6A,5") } ar®.[b, + (d&,/dd,)6r, |dé 


When these expressions are introduced in the right- 


hand side of Eq. (33) we obtain 


(r/a)(bo.’ — P,'P ) = 
(r/a)[@,’(d®,/dd,) — &,(d/dX,)(® VBA (36 


We now define the quantity @, = [d&,/d(\,£)] and we 


may write 


®,.’ = \,®, and d&,/dr\, = Ed, 


Introducing these expressions in Eq. (35) gives 
[r/a(o,’ — $)'6,) J = 
[r/a(tr.b,2 — bb, — £A,6,4,) ]B6A, (37) 


9 


Substituting Eqs. (35) and (37) in Eq. (33), 
through by 6\,, and taking the limit as 6A, > 0 gives 


dividing 
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é: Pret 
2; [ arb-dé = [(r/a) (EA. — Bb, — EX ®,,) |? 


Replacing #, and #, by their values now gives 


* 
2d | ar®dé = 
. 


$1 


te Ar a) (ED m2 — ©!’ — Em p ”) NE 
Finally, solving Eq. (26) for ®,” and substituting the 


result in the above expression gives 


N, = (1/2A,? 


JE (38) 


?.&,’ + aré&[A2 — (n?/r? 


The use of Eq. (38) will facilitate the calculation of 
the coefficients of series expansions in terms of the 
It should be noted that for given 
(3S 


functions &,(€). 
boundary conditions, the right-hand side of Eq. 
may be further simplified. 

We have that 
solution of Eq. (23) may be obtained in the separated 


now established the homogeneous 
and, because of the orthogonality of 


, the coeffi- 


form of Eq. (31), 
the characteristic functions ®,(£) and Z,(¢) 
cients A,,, may easily be obtained from the initial 
conditions. We now consider the nonhomogeneous 
solution corresponding to the function Q,(é, ¢, ¢) of 
Eq. (23). 


This solution may be taken of the form 


H,*é (60 = > > G, 


0é¢60 = % Dd K,,,()8(Z,(¢ 10) 
1 ¢ l 
Substituting Eqs. (39) and (40) in Eq. (23), and 


equating the coefficients of like terms in ?,7,, we 


obtain for G,,,, the differential equation 


G, mf + a(aA,;? + ¥ G, p —ak, pt 1] 


the particular integral of Eq. (41) gives the desired 


functions as 


Gis 


/ 


Thus, the complete solution of Eq. (23) may be written 


as 


(3) PARTICULAR CASES 


We have shown above that the determination of 
transient temperatures in shells of revolution may be 
reduced to the solution of Eq. (26) together with bound- 


ary conditions, such as those of Eq. (29). When the 
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SPACE 














set of functions ®,(£) has been determined, the solutions 
for H, is reduced to a routine matter. We now 
consider various particular cases and examine the 
resulting equation for ®,. 


Conical Shells 

It is convenient to take = y, where y is the distance 
from the apex of the cone measured along the middle 
surface (Fig. 2). It should be noted that at y = 0, 
the assumptions we made concerning the magnitude 
of ¢/r; and ¢/r, are violated, so that the present solu- 
tions would not apply; this has previously been noted 
by Hoff in connection with the theory of shells. With 
this choice for £, we have 


r=ysing, z= ycoss, a= 1 


where @ is the half-cone angle. Substituting these 


values in Eq. (26), we obtain 
b,” + (1/y)@,’ + [d,2 — (n?/y? sin? B)]®, = 0 (44) 


which is Bessel’s equation of order n/sin 8. The solu- 


tions are, therefore, of the form, 
©, = AJm(Asy) + BJ _m(Asy) (45) 


where m = n/sin 8. When m is a positive integer, 
J_m is not distinct from J,,; this is not important, 
however, since sin 6 is an integer only for 6 = O and 
8 = 2/2, which are not of interest here. Eq. (45) is 
not the complete solution for 7 = 0, which represents 
axially symmetric temperature and also enters into the 


case of even temperature distribution. For n = 0, 
the solutions ®, are given by 
&(y) = ASA) + BYol(Axy) (46) 
Spherical Shells 
In this case, we let £ = ¢ where ¢ is the inclination 


of the normal to the middle surface of the shell (Fig. 3). 
We then have 

r= Rsne, «= Rl —cos¢), a=R 
and Eq. (26) becomes 


®,” + (cos g/sin ¢)®,’ + 
[R2\,2 — (n2/sin? ¢)]®, = 0 (47) 


SCIENCES—DECEMBER, 1958 


It is convenient to change the variable to 7 = cos ¢: 


with this transformation Eq. (47) becomes 


(1 — n?)(d°@,/dn?) — 2n(d®,/dn) + 
{R2A,2 — [n?/(1 — 92)]}o, = 0 (48 


Eq. (48) is recognized as the associated Legendre 
equation with 


R>,? = s(s + 1) 


ig fT 2...) 


and the two basic solutions of the equations are the 
associated Legendre functions P,”(n) and Q,”"(n). The 
solutions for ®, are therefore in terms of the tessera] 
harmonics, and we have 


®,,(¢) =AP,"(cos ¢)+BQ,"(cos ¢) (49 
In particular, for axially symmetric temperature 
n = 0, and the functions reduce to the Legendre poly- 


nomials P,(cos ¢) and Q,(cos ¢). 

Similar analyses may be carried out for shells of 
other shapes. The governing equation for the func- 
tions ®,(£) is found to be related to well-known equa- 
tions of mathematical physics associated with various 
curvilinear coordinate systems. 


(4) CORRECTION FUNCTIONS FOR SLAB SOLUTIONS 


It is sometimes convenient to express solutions of 
Eq. (15) in terms of simple solutions of unidimensional 
heat-flow problems. The temperature in a shell of 
revolution is then expressed as a sum or product of a 
slab solution with a suitably defined correction function, 
which accounts for conduction in the meridional 
We will 
consider here two cases for which this further simpk- 


direction associated with the coordinate &. 


fication is possible. 
Case 1 
perature distribution with homogeneous-boundary con- 


Consider the case of axially symmetric tem- 


ditions of the type 


—k(OT/d¢)+A,T=0, (=;/2 
—k(OT/0¢)+A2T =0, t= 
—(k/a)(OT/OE)+Bi\T=0, §E=6, 
— (k/a)(OT/0t) + B.T =0, gé=€ 


It should be noted that Eqs. (50) allow for convection 
from the surface into a medium at zero temperature. 
Most convection problems can be transformed to this 
form by taking the environment temperature as a 
reference base, provided that this environment tem- 
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perature is uniform and equal for all faces. We assume 
further that the initial condition may be written in 


the form 
T(é, ¢, 0) =f(E) EO) (51) 
For this case Eg. (15) becomes 
1 /a2) {(0?7'/0E?) + [(r’/r) — (a’ /a) |(OT/dE)} + 
(0?7'/0¢?) =(1/a)(OT/Ot) (52) 
We now assume a solution of the form 
T(é, 6, D=W(E, t)W€, 0) (53) 
where II’(¢, f) is a slab solution satisfying the equation 


(0° W/0¢)?= (1 /a)(OW/ot) (54) 


together with the boundary conditions 


—k(00W/d¢)+4i1W=0, ¢=h/2 \ sien 
—k(OW/d¢)+A,W=0, ¢=—h/25 
and the initial condition, 
W(¢, 0) =2(¢) (96) 


Substituting Eq. (53) in Eq. (52) and making use of 
Eq. (54) reveals that the correction function W(&, ¢) 


must satisfy 


(] a”) }(O°'W £2) + [(r’/r) — (a’/a) |(OW dé)} = 
(1/a)\(OW/Ot) (57) 


The substitution of Eq. (53) in the boundary condi- 
tions (50) and the initial condition (56) shows that the 
form (53) satisfies all the conditions of the original 
boundary-value problem, provided that the correction 
function W(é, ¢) satisfies the following boundary and 


initial conditions 


—(k a)(OV 0£)+ By =0, &=£,| (58) 
P ID) 

—(k/a)(OW/0E) + BW=0, E=Ef 
w(é, 0) =/f(£) (59) 


We have therefore shown that, for the case of 
homogeneous boundary conditions, the temperature 
function 7(£, £, ¢) is the product of an appropriate 
slab solution with a correction function W(é, ¢), which 
must be a solution of Eq. (57), subject to the boundary 
conditions (58) and the initial condition (59). It 
should be noted that this function V may be interpreted 
physically as the one-dimensional temperature distribu- 
tion in a curved rod of shape defined by the parametric 
Eq. (1). This method is therefore an extension to 
curvilinear coordinates of the method mentioned by 
Carslaw and Jaeger.‘ 

Case 2—We now consider the case of axially sym- 
metric temperature when the shell is subjected to a 
time-varying heat flux on the surfaces (= +//2, the 
insulated. The 


boundaries =& and = being 


boundary conditions are 


—k(OT/OH)=qlt), ¢=h/2 } 
—k(OT/0f) =qglt), §¢=—h/2 (60) 
—(k/a)(OT/0E)=0, E=&, &o \ 


y=y, and y=, where y,+0. 


We further assume that the initial condition on 
temperatures may be expressed in the form 


T(é, €, 0) =f(E) + 2(0) (61) 
In this case, we investigate a solution of the form 
T(t, ¢, D=We, O+WEE, OD (62) 
where II’(¢, ¢) is a slab solution satisfying the equation 
0? /d¢? = (1/a)(OW/Odt) (63) 
subject to the boundary conditions 


OW jalt), ¢=h/2 | (64) 


— —_ 
or lq@(t), f= —A/2S 


and the initial condition 


We, 0) =2(¢) (65 


Substituting Eq. (62) in Eqs. (52), (60), (61), and 
making use of Eqs. (64) and (65), we now find that the 
function ¥(é, 4) must again be a solution of Eqs. (57 
In this case, however, the boundary conditions are 


—(k/a)(OW/0E)=0, E=&1, bo (66) 


and the initial condition is 


) (67) 


W(é, O) =f(E 


Thus, it is clear that, for many cases of interest, the 
temperature solutions may be obtained by a suitable 
combination of a slab solution and a correction function. 


(5) EXAMPLE 


As an example, we consider the transient temperature 
distribution in a conical shell subjected to time varying 
heat flux on its outer surface, all other boundaries 
being insulated. The initial temperature is assumed 
equal to zero. As a consequence of these conditions, 
the temperature distribution is axially symmetric and 


must, therefore, satisfy 


(027/dy?) + (1/y¥) (OT /Oy) + (027 /d¢2) = (1/a)(OT/2L) 
(68) 


where y is the distance along the cone middle surface 
from the apex. Since this shell theory does not apply 
at the apex, we consider the frustrum bounded by 
The boundary and 


initial conditions are 


k(0T/06) =q(y, D=PW)Qh, ¢=h/2 
(OT or) =0, ¢= —h 2 (69) 
(OT /Oy)=0, y=, Ve 

T(y, £, 0) =O (70) 


Note that the heat flux g(y,¢) is defined as positive 
when heat is transferred into the shell. 
The solution is assumed of the form 


T(y, §, D=Ay, &, O+C/2AR)Ch+S2) QO (y) (71) 


where // satisfies the conditions 
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OH/or=0, ¢=+h/2 } 
OH/oy=0, y=, y (72) 


Hy, ¢, 0)=— [2(0)/2hk ch te2yr(y) 9 


We may easily verify that Eq. (71) and (72) satisfy 
all the boundary and initial conditions on the tempera- 
ture. The function //(y, ¢, ¢) must also satisfy the 
differential equation 
(0°/T/Oy") +(1/y)(OlT/Oy) + (07/7/0¢?) = 
(1/a)(O/7/dt) — [Q(T (y) /hk]+ [(ch+6)2/2hk] X 

i (T'(y) /a](dQ/dt)+ Qir’+(1 y)I’}f (73) 

The homogeneous solution of Eq. (73) can easily 

be obtained by separation of variables; thus, we assume 


H = ®#(y)Z(¢) F(t) 


and find that #, Z, and F must satisfy Eqs. (44), 
In addition, ® and Z 
must satisfy the boundary conditions 


(27), and (28), respectively. 


&'(y)=0, v=, Ve 


Z'(¢)=0, ¢=+h/2 
The appropriate solutions for ®, are therefore 


we (74) 
, P ‘ 
&,=Ji(Axy) —R,Vo(dsy)  (s=1, 2...)9 
where the eigenvalues \, are the roots of 
Ji(As)¥1 Vi(Agy2) — Save) Vi(An) = 0 
and \y=0, Ry=Ji(Asvi)/ Vilas). The roots of the 
transcendental equation may be obtained graphically. 
The solutions for Z,(¢) are 
Zo=1 l 
Z (¢) =cos yp[(h/2)—¢], (p=1, 2...) 


(75) 


where y,»=(p2/h), p=1, 2... and yo=0. 
It is then clear from Eqs. (74)—(75) that the homo- 


geneous solution is of the form 


H(y, ¢, )=—[2(0)/2hk], 4 YO Ay, X 
s=0 p=0 
exp [—a(\.? +7 ,7)t]®.(y)Z,(¢) (76) 


where the coefficients A,, are obtained from the initial 
conditions (72); we have, therefore 


ye h/2 
f f (ChLEE) WT (y)O,(y)Z (Odedy 
4 = ¥1 h/2 
a ol ye feh/2 
ff yP,"(v)Z_7()dedy 
v1 h/2 


We now turn our attention to the nonhomogeneous 


(77) 


solution of Eq. (73). In this case, we assume a solu- 


tion of the form 
H*(y,60 = DY Dd Gyp(t)%.(y)Z,(0) 

s=0 p=0 

and we determine the functions G,,(¢) by substitution 

of this expression in Eq. (73). We must first expand 

the functions on the right-hand side of this equation 

into series in terms of the characteristic functions. In 


SPACE SCIENCES 


DECEMBER, 1958 


view of Eq. (77) we have 
[T'(y)/2hka ](th+¢2)(dQ/dt) = 


(1/2hka) 3) >> (d2/dt)A,,%,(y)Z,(¢ 
s 0 0 


D 
i 


[2(t)/2hk][T"+ (1/y)T’](Gh+e2) = —[Q(1)/2hk] x 


2. 2+. a hs OANZ AC 


ae 
0 p=0 


and we may also write 
IT(yQ)/hk] = [2 /hk]) SS DY Dy, O.0y)Z,(¢ 
s=0p=0 


where 


) 


vo feh/2 
[ | IT (9) ®()Z (OC) dedy 
~/ Vile h/2 _ 
V2 h/2 
f f y,"(y)Z,"(O)dedy 
My h/2 


Substituting all those expressions in Eq. (73) gives 
The 


Dy»= 


the differential equations for the functions G,,(f). 

result is 

Gig + a(rA,? + 5 Pg oe = (a hk)[D Tas 
(1/2)A,p,A.7JQ(t) — (1/2hk)A 


p2’(t) (79 


It should be noted that for the special case where 
s = p=0, Eq. (79) reduces to 


Goo’ = [aQ(t)/hk|Doo — (1/2hR)(dQ/dt)Amw (SO 


Since we have already satisfied the initial condition 
on //, we are only interested in the particular integrals 
for Eqs. (79) and (SO), and we therefore have 


t 
Goo(t) = (a hk) Dw [ O(r)dr — (1 Dhk)Apo, (Q(t) — Q(0 ] 
0 


G,,(t) = (a/hk)[D,» — (1/2)Aspd5?] X 


t 
few [— a(vA,2 + Yr lt — 7) JQ(t)dr — 
0 


} 


(1 2k)» exp [ — a(A.2 + yp?)(t — 7) JQ’(r)dr 
0 


These expressions may be combined into the single 
expression 


t 
Guplt) = J Faplt — 9) {(a/hb) Day - 
( 


) 
(1 2)A 5 pAs? ]Q( r) — (1/2hk)A ; pd" ( r){dr (S1 
where 
F,,(t) = exp[ — a(A.?+y7,7)], 
The complete solution for // is, therefore 
H(iy,60= 4% YX 4G,,(t) — [2(0)/2hk] x 
s=0p=0 


A .pF,,(t)} ®(y)Zp(¢) (83) 


Substituting Eq. (83) into Eq. (71) yields the tempera- 


ture solution. The second term of this expression may 


be combined into the series term by making use of 
Eq. (77); thus, we have 
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TRANSIENT 


Tyo. = LD De {G,,(t) — [2(0)/2hk] x 
5s = @ 


p=0 


A,,F. p(t) + 1/2hR)A QDI PZ (0) (84) 


It should be noted that, as pointed out in the previous 
section, this solution could also be obtained by the use 
of an appropriate slab solution and a correction func- 
tion. A numerical example is presented in Appendix B. 


(6) CONCLUSIONS 


We have shown that, by making suitable approxima- 
tions, which are based on geometric considerations, a 
simplified heat-conduction equation may be derived 
for shells of revolution. This equation has the ad- 
vantage of being separable for all shapes of shells, and 
leads to analytical solutions in terms of classical 
functions. In addition, for many cases of interest, 
the transient-temperature solution for a_ shell of 
revolution may be obtained in the form of a sum or 
product of a slab solution and a correction function. 
These techniques should prove useful in many problems 
which now necessitate the use of automatic computa- 


tion methods. 


APPENDIX (A) 


In the development of the simplified heat-conduction 
equation for shells, it was assumed that the term 
(l/r) +(1/r9) (OT 0) can be neglected in comparison 
to (0°7 O&). In order to check further into the 
validity of this approximation, we consider the tempera- 
ture distribution in a spherical shell and investigate 
the influence of the neglected term upon the solution. 

For a spherical shell subjected to uniform conditions 
over its surfaces, the temperature depends only on ¢, 


and Eq. (15) reduces to 
(0°77 /0¢?) + (2/R)(OT /0¢) = (1/a)(OT/Ot) (A-1) 
while Eq. (16) gives 
0?7/0¢? = (1/a)(OT/Ot) (A-2) 
We first consider the case of a shell subjected to a 


constant heat flux g over its outer surface, the inner 


surface being insulated. The boundary conditions 


are then 


OT /o¢ = (1/k)qg, § =h/2 


l - 
oT or = 0, c=—h 2{ (A-3) 


Using Eq. (A-1), the solution subject to the boundary 


conditions (A-3) is easily obtained as 


T(¢, t) = (q/2hk) [fh + ¢?] + > ALA 
>=0 
(A-4) 


where the functions Z, are solutions of the equation 


Z," + (2/R)Z,' + 7p°Zp = 0 (A-5) 
with Z,'(+h/2) =0 


We find, therefore, 
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Z=1 ] 
Z,=0e~* ® [cos (pre h) + L, sin (prt /h)], p=1,2...4 
A-6 
where L,»= —prRh, podd 
=h/(prR), peven 


The exponential time factor 7,° is given by 


Vp? = (p?x?/h*) + (1/R? (A-7) 


Using Eq. (A-2), the solution of the corresponding 
problem is easily found to be 


T(¢, t) = (g/2hk)(¢h+¢7) + >. Awe dye A-S 
where 
Zo=1 | 
' > (A-9 
Z, = cos (par/h)[(h,/2)-—¢], p=1,2...5 
and Vp” = (p?x?/h? (A-10) 


It is observed that the first term in Eqs. (A-4) and 
(A-S) are identical. In addition, comparison between 
Eqs. (A-7) and (A-10) shows the percentage error in 


the time exponential factor to be 


(Ayp?/ Vp?) = (h/R)7(1/p?x?)/[1 + (8/R)2(1/p? x?) | 
(A-11) 


These errors are, obviously, very small for the values 
of h/R appropriate to thin shells. Finally, by using 
Eqs. (A-6) and (A-9) we may compute corresponding 
mode shapes for Eqs. (A-1) and (A-2), it is then found 
that the temperature variations with ¢ predicted by 
Eq. (A-2) differ negligibly from those predicted by 
Eq. (A-1). It is, therefore, concluded that our as- 
sumption is perfectly justified for this type of boundary 
conditions. 

We now consider the case of a spherical shell with 
temperature 7’, prescribed on the outer surface, the 
inner surface being maintained at zero temperature. 
The boundary conditions are then 


T[(h/2),t] =TA 


_19) 
T[ — (h/2), t] =0 $ alate 


Using Eq. (A-1), the solution is found to be 


T(¢,t) = (T,/2h)(h + 2¢) + ¥ ApZp(b)e7%" (4-13) 
l 


hb = 


with Z, = @ */® cos (2p + 1) (me /h) (A-14) 


and Yo? = (2p + 1)?(w?/h?) + (1/R?) (A-15) 


ip 


From Eq. (A-2) we obtain the solution 


T(k, t) = (T/2h)(h + 26) + DS A,Z,(He~" ~(A-16) 
t l 

where Zp»=cos (2p+1)(re/h) (A-17) 

Y p= (264+1)7(9?/h?) (A-18) 


Comparison of Eqs. (A-13)—(A-15) with Eqs. (A-16) 
(A-18), respectively, reveals that, just as before, 
there is a negligible difference between the two solutions 
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Fic. 4. Comparison between slab solution and present solution 


for the values of i, R appropriate to a thin shell. We 
conclude, therefore, that the assumption made is justi- 
fied for a wide range of boundary conditions. 


APPENDIX (B) 


This numerical calculation is included to indicate 
the difference between results of the present method 
and the slab solutions. 

The example selected is an aluminum conical shell, 
as shown in Fig. (2) with y,=1 ft., ye=5 ft., and h= 
0.025 in. 
solution is given by Eq. (84), where %,(y) is defined 


For axially symmetric temperature, the 


in Eq. (74). The first four eigenvalues for this case are 
A,=0, 0.842, 1.615, 2.380, for s=0, 1, 2, and 
Z,(¢) is given by Eq. (75) with y,=0, 1.51, 3.02, 4.5 
6.04 for p=0, 1, 2,3. Similarly F,,(¢) is obtained from 
Eq. (30). 

In this case, we assume a heat flux to the surface 


9 
o. 
9 
0, 


{= 4/2 as 


q(y, t) = Q(T (y) = QU — ky*)e-4 (B-1) 


where Q, x, and ¢ are constant. 
A,, and D,, are obtained by substituting Eq. (B-1) 

into Eqs. (77) and (78). 

A,»=K,E, for s=0, 1, 2... and p=0, 1, 2... 


where Ko= 1— (x/2) (91? +92?) (B-3) 


f y( 1 — xy?) [Jo(Asy) — Rs Yo(Asy) dy 
Kw 2h. 


f y[Jo(Asy) — Rs Yo(Asy) Pedy 


¥1 
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de [(Asy)? }Jo(As¥) — Rs Vo(Asy) } 3? 
Ast [vy {Jo(Awy) — R, Voy) } 223? 


K.= 


M/s 
} (ch + £7) cos yp[(h/2) — ¢]dé 
ae 


. ) } 
E, as. 8/2 7 = 
| cos” ¥Y rl(h 2) — ¢|d¢ = 
b=1,2 B-6 
Di, =K,M,, s= 0,1, 2 , ~=O,1,2 B-7 
where Mo=1 B-8 
cos ¥ pl(h 2)-—¢ ld¢ 
vw, = — =U p - 
cos? y,[(h/2) — ¢]d¢ 
B-9 


Finally, for the example, Eq. (S4) becomes 


Tiy,6,0= > DY §S[Oa/hk]/ [a(a2+y,2)—e]! Xx 
) 


0 p { 
[exp (— et) — exp {—a(A,?+y,2)t} ]X 
[Dept DA sp¥p7]Os(W)Zp(O) — (B-10 
For the numerical calculation, we take 


= ().03 ft.-? 


K 

e = 0.0167 sec. 

t = 60 sec. 

a = 9.7X 10-4 ft.?/sec. 


hk = 7.13X10-* B.t.u./sec.-°F. 


The results of the calculation are presented in Fig 
t, where the normalized temperature of the shell 
at this value of time the inner and outer surfaces are 
at the same temperature for plotting purposes—are 
plotted as a function of y. The dotted line is the cor- 
responding temperature distribution obtained from the 
slab solution where the heat flux at any y corresponds 
to that given in Eq. (B-1). The value of temperature 
for y=, and ¢= —//2, calculated by the slab solution 
method, is taken as the norm. It is observed that 
large discrepancies can occur when meridional conduc- 
tion is neglected, the slab solution giving a temperature 
24 per cent too low at y= yp» in this example. 
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Turbulent Heat Transter Through a Highly 
Cooled, Partially Dissociated Boundary Laver' 


Pretrer H. Rose,* RONALD F. PROBSTEIN,** AND Mac C. Apams*** 


y 1VCO Research Laboratory 


SUMMARY 


The problem of heat transfer from high-temperature air 
through turbulent boundary layer to a cold surface is con 
sidered both analytically and experimentally Heat-transfer 
data obtained in shock tubes are presented and correlated by a 
semiempirical theory which includes the effect of atomic diffusion 

The distinguishing characteristics of turbulent boundary 
ive with dissociation and large cooling are considered It is 
shown that the equations governing such flow, after certain 


roXimMations, can be represented in a form similar to the classi 
cal equations for a turbulent boundary layer 

\n approximate theory is proposed for turbulent heat transfer 
fora highly cooled boundary layer on portions of the body where 
ssure gradient is negligible in the case of blunted bodies 


th pr 


the pre 
of revolution in high-speed flight 

Experimental results obtained on the cylindrical portion of a 
hemisphere-cvlinder model are presented for conditions simulat- 
sec., where up to 30 per cent of 
the molecules are dissociated. Reynolds Numbers of 2.5 & 10°, 
based on local fluid properties external to the boundary layer, 
wert Reynolds Number and 


flight speed were not obtained simultaneously, due to structural 


ing flight speeds to 21,350 ft 


ichieved The larger values of 
shock tubes; 
i way that the important effects of each could 


limitations of the however, the experiments were 
conducted 1n such 
be determined 

In the experiments the Mach Number external to the bound 
iry laver varied between 1.7 and 2.2. The corresponding Mach 
Number for blunted nonslender bodies in flight would have a 


maximum value between 2.5 and 4; however, it is shown that 


these differences in Mach Number are not important for such 
bodies 
SYMBOLS 
( = skin-friction coefficient, 7/(1/2)p,.u.? 
( = incompressible skin-friction coefficient 
‘ = mass fraction of ith component 
é = specific heat per unit mass at constant pressure for a 
mixture 
C; = specific heat per unit mass at constant pressure of 
ith component 
D = laminar diffusion coefficient 
DT = thermal diffusion coefficient 
Dr = turbulent diffusion coefficient 
= internal energy per unit mass 
h = atomic dissociation energy times atom mass fraction 
} = enthalpy per unit mass of the mixture, including 
dissociation energy 
h = perfect gas enthalpy per unit mass of ith component 
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h = heat evolved im the formation of component at 
O°K., per unit mass 
h = perfect gas enthalpy per unit mass 
H = total enthalpy, h + (u?/2 
k = thermal conductivity 
L = Lewis Number, Dp@,/k = Dpo u 
Lr = turbulent Lewis Number, Drp@,/« 
MJ. = Mach Number external to boundary layer 
Vf. = Mach Number of moving shock in shock tube referred 
to speed of sound in quiescent gas in front of it 
Vu = Nusselt Number, gxvo/u,( 7, h 
f: = pressure 
t = initial pressure in shock tube 
) = surface heat-transfer rate 
evlindrical radius of body 
R = gas constant of 7th component 
Re = Reynolds Number, p,4.x/u 
St = Stanton Number, g/p,.u,(H h 
1 = absolute temperature 
= x component of velocity 
= ycomponent of velocity 
= diffusional velocity of 7th component 
= mass rate of formation of component per unit 
volume and time 
\ = distance along meridian profile of body 
; = distance normal to body surface 
3 = exponent in Lewis Number correction factor, Eq. (82 
4 = average specific heat ratio 
€ = eddy viscosity, Eq. (8 
A = eddy thermal conductivity, Eq. (9 
V1 = absolute Viscosity 
v = kinematic viscosity 
p = mass density 
0 = Prandtl Number 


= surface shear 


Subscripts 


é = external flow conditions 

1 = ith component of mixture 

stag = stagnation point conditions 

w = wall conditions 

x, ¥ = differentiation with respect to indicated variable 
Superscripts 


= time-averaged quantity 
= fluctuating quantity 
, = conditions corresponding to reference enthalpy 


(1) INTRODUCTION 


AN HYPERSONIC SPEEDS, as in subsonic and super- 
sonic flow, if the Reynolds Number is high 
enough the boundary layer becomes fully turbulent 

that is, the local flow is no longer steady and the 
velocity components fluctuate in a random manner 
Since the resultant turbulent will 
cause an increase in shear, heat transfer, and diffusion, 


mixing process 


it becomes of great importance to determine these 
quantities for the case of atmospheric hypersonic 
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flight, as for example in the re-entry of a ballistic missile. 

Because of the attendant high heat-transfer rates 
associated with atmospheric hypersonic flight, all 
bodies must be blunted to some extent. Furthermore, 
because of material limitations, there must be a large 
degree of wall cooling present under hypersonic flight 
conditions, since the wall temperature is necessarily 
small (generally below the dissociation limit for air) 
in comparison with the temperature of the air which 
has passed through a detached normal shock wave at 
high relative velocities.| Thus, for example, for 
ballistic missile re-entry conditions, the temperature 
in the shock layer may be as high as 8,000°K. Be- 
cause of the conversion of the large amount of kinetic 
energy associated with the high flight velocities into 
thermal and chemical energy, the ratio of the wall 
enthalpy (,,) to the free-stream stagnation enthalpy 
(/7,), and the ratio of the external density (p,) to the 
wall density (p,) are small in comparison to unity. 
The density ratio p,’p, is small not only because 
Felts 
the air external to the boundary laver is less than that 
of the air at the surface. Furthermore, since the 
actual Mach Number cf the flow external to the 
boundary layer for a nonslender blunt body will in 


<_ 1, but also because the molecular weight of 


general not exceed 3 or 4, the Mach Number effects, 
usually associated with hypersonic flows, will not be 
large, although the air temperature involved in the 
boundary layer will be high. 

Of course, because of the extremely high tempera- 
ture in the boundary layer the air will be dis- 
sociated, and it becomes necessary to consider the 
diffusion of atoms through the boundary layer in any 
calculation of the surface heat transfer. In order to 
simplify the calculations, in the present paper we 
will assume that the flow in the boundary layer is in 
thermodynamic equilibrium. That is to say, we will 
assume that the reaction times are short in comparison 
with the time of convection or diffusion through a 
boundary layer, so that chemical equilibrium prevails 
throughout. Thus, either the temperature or con- 
centration distribution along with the pressure is a 
sufficient description of the thermodynamic state of 
the boundary layer. For laminar boundary layers 
at high densities, this has been shown to be an excellent 
approximation (e. g., see references 2 and 11). 

Although recently 
achieved”: *: !* in calculating the laminar heat-transfer 


some success has been 
rates on highly cooled blunt bodies for hypersonic flight 
conditions, no such straightforward analytical approach 
can be applied to the treatment of the turbulent prob- 
known, the Prandtl 


boundary-layer concepts are still applicable to a 


lem. However, as is. well 
turbulent flow, and a mathematical treatment can, in 


+ In this report the primary interest is in bodies sufficiently 
blunted such that the air which enters the boundary layer has 
passed through a nearly normal shock wave. Thus, the inviscid 
flow external to the boundary layer has approximately constant 
entropy and is therefore nearly isentropic. Such bodies are termed 


nonslender blunted bodies. 
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SCIENCES 1958 
principle, be carried out by assuming the motion to 
be separable into a mean flow and a superposed turbulent 
fluctuation, the mean value of which is zero. On 
this basis one can, at least formally, derive the appro- 
priate differential equations of mean motion, energy 
and continuity for the turbulent boundary layer in aq 
compressible dissociating flow. Apart from the at- 
tainment of general flow characteristics and certain 
specific results from these equations, it is, in general 
not possible to develop any exact analytic scheme for 
the ccmputation of the heat-transfer rate. Therefore, it 
is necessary to rely upon semiempirical theories and 
experiment in order to obtain the desired results. 
In the following sections we will attempt to outline 
a general theoretical, empirical, and experimental 
attack for determining the turbulent heat transfer 
through a highly cooled, partially dissociated boundary 
layer 


(2) THe TuRBULENT BOUNDARY LAYER IN A 
DISSOCIATING GAS 


(a) Equations of Mean Motion 


The equations for the mean motion of a dissociating 
turbulent boundary layer can be derived from the 
general equations of motion, energy, ete. These 
general equations may be found in references | and 2 
and are reproduced here, after the usual boundary- 
layer simplifications for the axisymmetric case, as 
follows. 


The continuity equation for any species: 


[O(pc;) /Ot] + (1/7r)(0 Ox) (purc;) + 
(1/r)(O0/Oy) (prrc;) = (0 Oy) [Di p(Oc; Ov) 4 
(D;"pc;/T)(OT /Oy)] + w; (1) 


The summation of Eq. (1) over all species leads to the 
usual form of the continuity equation 


(Op Ot) + (1 r)[O(pur) Ox] + 
(1/r)[O(pzr) Ov] = 0 (2 


The equation of mementum: 


[O(pu) Ot] + (1/r)[O(pu?r) Ox] + (1 7)(0 Oy) X 


(puvr) —(Op Ox) + (Ody) [u(Oudy)| (3 


The energy equation, neglecting thermal diffusion,{ 


becomes 


[O(p/Z) Ot] + (1/7r)(0/0x) (pur) 
(1 ‘r)(0/Oy) (pvlTr) (0 Oy) [(k &,)(Ol/Oy) | + 
z, [Dip — (k Cp) (hi — h;®)(Oc;/Ov) + 


Ll 


(0/Ov); (u/2)[1 — (1/e)](Ou? /dy)} (4 


It should be noted that the enthalpy /7 includes the 
energy of dissociation and, in addition, the kinetic 
energy- thus, 

tIn the case of a boundary layer in thermodynamic equilib 


rium, it can be shown that thermal diffusion is unimportant 
for temperatures less than 10,000°K. 
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TURBULENT 


H Veil + RT —h®) + (uv? 2) 
Dei(h — ho) + (y2/2) 
ho thp + (u?/2) =h 4 (u2/2 
where 


h Yc (e, + RiT) = doch 
: 


i 
and the dissociation energy is given byt 


hy — Sieh ’ 


The boundary layer is assumed to be in thermo- 
dynamic equilibrium, and this, in effect, furnishes the 
“equation” of state. 

In the analysis to follow, both the Prandtl Number 
é,u/k and the Lewis Number Djpé,/k will be taken 
as unity. Semiempirical corrections will be made 
later to account for values different from unity. In 
the case of a laminar boundary layer, the Lewis Num- 
ber, which is the ratio of the particle to thermal diffu- 
sivity, is believed to be approximately 1.4 in the 
temperature and pressure ranges being considered 
(e.g, see reference 2). It might be expected that 
an effect of the laminar Lewis Number will be 
present for the turbulent boundary layer, since, as 
will be shown later, the free-stream atom concentration 
may persist well into the boundary layer and even to 
the extremities of the laminar sublayer. Thus the 
diffusion of atoms through the laminar sublayer would 
be similar to the diffusion in the case of a purely 
laminar boundary layer. It will be seen in Section (3) 
that the heat-transfer measurements are best corre- 
lated with a laminar Lewis Number correction factor. 

No information exists for the magnitude of the 
turbulent Lewis Number; however, its value is prob- 
ably nearer to unity than is the laminar Lewis Number, 
as in the case of the Prandtl Number. 

Following the usual procedures,’ the fluid properties 
are written as a sum of a time-averaged quantity plus 


a fluctuating quantity—thus, 


uw=iuitu, h=h+h’", p=pt pp’ 
and the time-averaged differential equations for the 


steady-state turbulent boundary layer become 
Continuity 


(0/Ox) (pur) + (0/Oy) (pir) + (0/Oyv)(p'v’r) = 0 (5) 
Momentum 


pu(Ou/Ox) + (pd + p’v’)(Ou/Oy) = 
—(Op/Oxv) + (0 Oy) [u(Om Oy) + €(OH Oy)| (6) 


Energy 
pu(o0H/dox) + (po + pv’) (0 ov) = 
(0/dy) [(k é,) (OH Oy) + (x /é,) (O/T Oy)] (7) 
The eddy viscosity ¢ is defined in the usual sense 


t If the reference level for enthalpy is chosen to be the ground 
state of the molecules at zero °K., then /;® for the atoms will 
be minus the dissociation energy and zero for the molecules. 
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e = —[pu'v’ (On Oy)] (8) 


Some explanation is required for the eddy thermal 
conductivity «x. In Eq. (7), x has been taken to be 
k = —[po'h’/(1/é,)(Oh dy) ] (9) 
This definition of the eddy thermal conductivity is 
consistent with its usual definition in the case of a gas 
with constant specific heat and no dissociation. Now, 
in the case of a dissociated boundary layer, the energy 
transported by the diffusion of atoms must be ac- 
counted for. 
The energy transported by the ordinary laminar 


conduction and diffusion processes is given by 
Energy transported across a plane normal to the 


main flow = k(OT dy) — dopvici(h; — h,®) (10) 


t 
where 7; is the diffusional velocity 
v, = —(D;‘c;)(Oc; Oy) —(k/ pci€p)(Oc;/Ov (11) 
for Lewis Number of unity. 
Eq. (10) then becomes 
Energy transported 
(k €,)>, [e(dh,/dT)(OT Oy) + (hy — h,°)(Oc; Oy)] 
1 
= (k c,)(Oh/Oy) (12) 
where use has been made of the relation 


Cp = Yc (dh; dT 


Thus, the energy transported by conduction plus 
diffusion is included in the term (k‘c,)(Oh Oy) for the 
case of a laminar Lewis Number of one. 

Now, in the case of turbulent fluctuations, the 
energy transported across a plane normal to the main 
flow is given by 


— ph’ = —pl>-é, vk; + Yoo" "th, — h&)] 
t 


According to the mixing-length concept 


Leah? = — Dede’ oé,)(h,/y) 


and 

Yv'c/(h; — h) = —YDz(0é;/dy)(h; — h 

1 1 
where D,; represents the turbulent diffusion coefficient. 
Thus, it follows that 
—po'h’ = (x/é,)(0/dy)[De(h, — hi®)] + 

1 
(x €,)(Lr — 1) >> (dé, oy) (h, — h;°) 


Ly, denotes the turbulent Lewis Number 


Now, with the approximation of the turbu- 


where 
Dypep/k. 
lent Lewis Number equal to unity, 
— po'h’ (x €,)(0 Oy) Deh; — he) = 

: 

(x €,)(Oh/Oy) (13) 


The importance of density fluctuations is not known. 
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Fluid property distributions through a highly cooled 
turbulent boundary layer. 


Some attention has been given to their influence 
(e.g., see reference 9), however, in the remainder of 
this report they will be neglected, although there is 
no clear justification for doing so. 

If now the density fluctuations are neglected, Eqs. 
(5), (6), and (7) become 


(0 Ox) (pur) + (0/dy) (par) = O (14) 


piu(Onu/Ox) + po(Ou/Ov) = —(Op/dx) + 
(0/Ov) [(u + €)(O%/Oy)] (15) 


pu(Oll/ox) + pa(d/7/dy) = 
(0/Oy)[(u + €)(O0H/dy)] (16) 


where the laminar and turbulent Prandtl Numbers 
have been taken as unity, 


Cpu/k = Cpe/x = 1 


(b) The Crocco Integral and Reynolds Analogy 


From Eqs. (15) and (16) it is apparent that the 
familiar relation between energy and velocity (the 
Crocco integral) still exists in the case of a dissociated 
boundary layer with zero pressure gradient, and with 
both laminar and turbulent Prandtl Numbers of unity, 
thus,t 


H = hy + (He — hy) /ue|u (17) 


where the subscripts w and e denote, respectively, the 
conditions at the wall and at the outer edge of the 
boundary layer. Here, the total enthalpy //, is 
equal to the stagnation enthalpy. 

It is of interest to examine the effect of wall cooling 
on the distribution of fluid properties through the 
boundary layer. Considering the case of zero pressure 
gradient and Prandtl Number one, the gradient of 
enthalpy through the boundary layer follows from 


Eq. (17) and is 
oh/dy = }((H, — hy)/u.] — ut (Ou/dy) (18) 


Now for most cases of practical interest the velocity 


gradient is positive, consequently, the enthalpy 
gradient is positive if 
(7, — hy)/u.> u (19) 


+ The bars (~) will be dropped hereafter. 
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and this condition is satisfied throughout the boundary 
layer, providing 


(H, — hy)/ue> u, 20) 


Rearranging Eq. (20), the condition for the maximum 
enthalpy to occur at the outer edge of the boundary 
layer is then 


.?/f1, < 1 = e/f,) 2] 


provided, of course, h,,/h, <1. In general, other fluid 
properties such as density, temperature, etc., will follow 
the same trend as the enthalpy, and they too will have 
no maximum point, providing Eq. (21) is satisfied. 

For a blunted body (such as a hemisphere-cylinder 
in high-speed flight, the ratio of the kinetic energy of 
the inviscid flow to stagnation enthalpy is never 
greater than approximately one half, thus, u,?/ 77, < 1, 
Then, in the case of a highly cooled boundary layer 
where the ratio of wall enthalpy to stagnation enthalpy 
is negligibly small, maximum values of the fluid proper- 
ties will occur at the boundary-layer extremities. 

If Eq. (21) is roughly interpreted in terms of a 
mean ratio of specific heats 7, then, in more familiar 
terms, 


,* 2(7¥ — 1)M,? h,, 
i, 2+ (7 -D)M-2 


The condition u,?/H, = 1 leads to M,? = 2/(¥7 — 1), 
which gives, for y = 1.4, M, = 2.24, and for 7 = 1.2, 
M, = 3.16. The latter value of 7 is more realistic for 
flight Mach Numbers above approximately 10. 

Some typical distributions of enthalpy, temperature, 
density, and dissociation are shown in Fig. 1 for the 
case of a turbulent boundary layer with an assumed 
velocity profile of the form u/u, = (y/6)'". The 
nature of the distributions is not critical to the 1/7 
power law, and they retain the same general character 
provided the velocity profile is of the turbulent type— 
i.e., a full profile. 

Since the fluid property distributions are not appre- 
ciably affected by Mach Number up tc a value of 
3.16 and since for a nonslender blunted body in flight 
the maximum value is approximately 4,{ then the 
Mach Number as such is not expected to be an im- 
portant parameter for turbulent heat transfer. The 
decreasing importance of Mach Number with increased 
cooling is also noted by Sommer and Short.'” 

Reynolds analogy between heat transfer and skin 
friction follows directly from the Crocco integral Eq. 
(17) or (18)—+.e., 


G/ pete. — hy) = T/ Pele? = Cy/2 (23) 


The left-hand side of Eq. (23) will be defined as the 
Stanton Number—i.e., 


St = g/pul(H. — hy) = ¢;/2 (24) 

t This value of Mach Number is obtained by considering an 

isentropic expansion back to the free-stream pressure for the 
flow that has passed through a normal shock. 
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It should be pointed out that Stanton Number is also 
commonly defined with recovery enthalpy H/, in place 
of the stagnation enthalpy /7, in Eq. (24). The well- 
known relation between Stanton Number, Nusselt 
Number, and Reynolds Number will also be used in 
future sections. 


Nu = St Rezo = qxo/u,(H, — h, (25) 

For the case of the Prandtl Number different from 
unity, the ratio of the Stanton Number to one-half 
the skin-friction coefficient has been found by previous 
authors to be a factor greater than unity. Colburn‘ 
determined experimentally for low velocities and small 
temperature differences that this factor was apparently 
the same as in the case of a laminar boundary layer— 


ie.. a laminar Prandtl Number correction factor 
correlated his data. Thus the relation Colburn 
proposed was 

St (Cy 2) gq — 2/8 (26) 
and with o = 0.72 the factor o—®/*®) becomes 1.25. 


Rubesin,® using an analytical approach, arrived at a 
factor which also depends on the laminar Prandtl 
Number but varied between a value of 1.17—1.20. 
Since the differences between Colburn’s and Rubesin’s 
results are small, the former result, being simpler, will 
be used in the present report to correlate the experi- 
mental data. Eq. (26) will be referred to as the 
modified Reynolds analogy, which has common usage 
in the literature. We also note that Eq. (26) is often 
written in the form 


St = (tc; Cri) (Cy 2) q— (2/8 (27) 


where cy, is the incompressible skin-friction coefficient. 
Although there is some arbitrariness in the selection of 
wall or external conditions for the evaluation of the 
Prandtl Number, the best available estimates indicate 
a small variation with temperature up to 8,000°K. 
Hence, the choice of temperature for its evaluation 
makes no significant difference. 

The nearly constant density distribution through 
the highly cooled boundary layer (see Fig. 1) suggests 
that the skin-friction and heat-transfer mechanisms 
are similar to the case of an incompressible turbulent 
boundary layer. While no rigorous argument can 
be made that the incompressible coefficient should 
apply, a plausible empirical argument can be made. 
The large changes in density and other fluid properties 
occur in the laminar sublayer. Further, since the mag- 
nitude of the momentum in the sublayer is a small frac- 
tion of the total boundary-layer momentum, as in the 
ordinary incompressible case, the overall growth of the 
should not be significantly 
Thus, it is reasoned 


thickness 
influenced by the sublayer. 
that the incompressible coefficients should apply for 
both skin friction and heat transfer. It will be seen 
that this assumption leads to a good correlation of the 


momentum 


heat-transfer data. 
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(c) A Correction for the Laminar Lewis Number 


For the case of a dissociating turbulent boundary 
layer, we have noted (see Fig. 1) that appreciable 
atom concentrations appear to persist to the extremities 
of the laminar sublayer. It may, therefore, be 
expected that the laminar heat-transfer mechanism 
through the sublayer will be similar to the case of a 
dissociated laminar boundary layer. Thus, analogous 
to the Prandtl Number correction, a laminar Lewis 
Number correction, for values different from unity, 
is proposed for the case of a dissociated turbulent 
boundary layer. The form of this correction can be 
most easily seen by considering the heat flux crossing 


any plane in a laminar boundary layer, 
gq = k(dT/dy) + Daph*(dcs dy) (28) 


Here h° is the dissociation energy per unit mass of 
atomic products, D4 the atomic diffusion coefficient, 
and c, the atom mass fraction. In deriving Eq. (28) 
the specific heats of the atoms and molecules have 
been assumed equal, which is approximately true for 
air in the high-temperature region in which the atom 
concentration becomes significant. If we express the 
heat flux in terms of the enthalpy we may rewrite 


Eq. (28) as 


ae (kb é,)(L _— 1) x 
[d(h°c4)/dy| (29) 


q = (k Ep) [d(ho + h%c,) dy 





The first term is the heat flux determined by the 
chemical enthalpy, which is the perfect gas enthalpy 
plus the enthalpy of formation. We then calculate 
the ratio of the heat-transfer rate for Lewis Number not 
equal to unity to the heat transfer with L = 1, as 


d(h°c4)/dy e 
) . 
d(hyg + h°c,4)/dy 


dis 


qr=1 


=~1+(L-1 


Taking the simplified one-dimensional case of no 
convection, in which we have linear distributions 
between a high-temperature surface and a cold wall, 


Eq. (30) gives 


1+ (L — 1)(hp,/H,) (31) 


Gus Qu=1 


where /ip, is the energy in dissociation at the high 
temperature surface and //, is the total enthalpy at 
this surface. 

On the other hand, for highly cooled 
boundary layers, both Fay and Riddell’? and Lees® 
have shown, for the case in which a large fraction of the 


laminar 


energy is carried by dissociated atoms, that the Lewis 
Number correction enters in the form 


+ (L? — 1)(hp,/H,) (32 


Givi Jr=1 = | 


The numerical calculations in reference 2 show that, 
at a stagnation point for an equilibrium flow, 6 = 0.52, 
and, for frozen dissociation, 8 0.63. It should be 
pointed out that the results of these numerical cal- 
culations, using these factors, were in good agreement 
Based on more 


approximate calculations, Lees® estimated that 8 =2/3 


with experimental measurements."! 
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for the frozen highly cooled flat-plate case. In all 
of these cases, the change in the exponent from unity 
results from the convection present in the boundary 
layer. 

Now a correction factor, depending on the laminar 
Lewis Number, of the same form as Eq. (32) may be 
applied to the turbulent dissociating boundary-layer 
heat-transfer problem. Of course, the determination 
of the exact value of 8 requires recourse to experiment, 
and in Section (3) the value obtained from the present 
experiments will be discussed in connection with the 
correlation of the data. 

Summarizing the results of Eqs. (25), (26), and (32), 
the modified Reynolds analogy may be expressed in 
terms of a Nusselt Number as 


— hy) = (c,/2)[1 + (L8 — 1) X 
(hp,/H,)]o./*Re, (33) 


Nu = qax pu, (LT, 


Based upon the previous empirical arguments. the skin- 
friction coefficient is taken equal to its incompres- 
sible value; thus, 


c,/2 = 0.029 (Rez)? (34) 


The Lewis Number correction term [in brackets in 
Eq. (33)] may now be interpreted as the ratio of the 
compressible to incompressible skin-friction coefficients; 
thus, 


Cy/Cy = 1+ (L° — 1)(hp, HT,) (35) 


This interpretation will be utilized in the correlation 
of experimental data. 


(3) TURBULENT HEAT-TRANSFER MEASUREMENTS FOR 
A HIGHLY COOLED DISSOCIATED BOUNDARY LAYER 
The shock tube and experimental techniques used 

in the present experiments have been described in a 

previous paper!! concerned with heat-transfer experi- 


SCIENCES 


DECEMBER, 1958 


ments under stagnation conditions simulating hyper. 
sonic flight. 

The stagnation conditions attained in the shock 
tube operation and their corresponding simulated 
Complete 
simulation of specific flight conditions over a given 


flight conditions are given in Fig. 2. 


geometry cannot be attained without simulation of 
all the parameters—Reynolds Number, Mach Number 
and enthalpy. As a result, in making shock tube 
experiments away from the stagnation point, a some- 
what different simulation philosophy must be adopted 
than was used at the stagnation point.” 

Shock tube experiments can be performed in regions 
where the properties of the flow are very closely 
matched to those actually encountered in hypersonic 
flight for blunted nonslender bodies. The Reynolds 
Number, the enthalpy ratio across the boundary 
layer, the degree of dissociation of the air at the edge 
of the boundary layer, and the Mach Number of the 
flow at the edge of the boundary layer, attainable in 
a straight shock tube, cover, for example, the ranges 
of these same variables encountered by ballistic missiles 
during re-entry remarkably well. 

A map indicating the region of Reynolds Number 
vs. free-stream-to-wall enthalpy ratio covered by the 
present experiments is shown in Fig. 3. 
constant temperature in the external flow (and tem- 


Lines of 


perature ratio across the boundary layer) are indicated 
on this Figure. In the region below the point at which 
natural transition occurs, a trip wire was used to 
produce turbulent flow. 

The Mach Number at the outer edge of the bound- 
ary layer varied between 1.7 and 2.2 in these experi- 
ments, while in flight this Mach Number would be 
less than 4.0, as discussed in Section (2). Though the 
present experiments do not include a large Mach 
Number variation, it is felt that this is not important 
in view of the discussion given in Section (2). 

The degree of dissociation under which these turbu- 
lent heat-transfer measurements have been made is 
shown in Fig. 4, plotted against Reynolds Number. 
Lines of constant flight velocity for equivalent stagna- 
tion enthalpy are also indicated. 

Having shown that the shock tube experiments 
cover a range of conditions broad enough to encompass 
the variations encountered in re-entry flight where 
dissociation is important, it becomes necessary to be 
able to correlate the results of these experiments with 
our previous theoretical concepts in order to predict 
This simulation philosophy 
4.e., that cor- 


turbulent heat transfer. 
has been adopted in the present program 
relation of experimental results, obtained under condi- 
tions similar to but not exactly the same as those en- 
countered in flight, allows the turbulent heat-transfer 
rates in flight to be predicted. 

The experimental turbulent heat-transfer program 
has been divided into four phases. First, it was 
necessary to establish the fact that both laminar and 
turbulent heat-transfer measurements could be made 
in a straight shock tube and transition could be 
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achieved. Second, the Reynolds Number range ac- 
cessible in the shock tube was covered to establish 
general turbulent heat-transfer trends. Third, ex- 
periments were performed to establish the effect of 
the variation of the fluid properties, including dissocia- 
tion effects, through the boundary layer. Fourth, 
an attempt was made to establish the effect of pressure 
eradients and body shapes by making heat-transfer 
rate measurements on various configurations. The 
fourth phase will not be covered in the present paper. 


(a) Laminar and Turbulent Heat-Transfer 

Veasurements 

In the first phase of this investigation, heat-transfer 
rate experiments were performed over a wide range of 
Revnolds Numbers, approximately 10' to 3 X 10°, 
on the aft portion of a hemisphere-cylinder body. 
The pressure distribution on this geometry was pre- 
viously measured employing the techniques, described 
in reference 11, of measuring the Mach angles created 
by weak disturbances on a body (see Fig. 5). These 
measurements indicate that the pressure on the cylin- 
drical portion is approximately equal to the free-stream 
static pressure. Furthermore, the Newtonian theory 
without centrifugal correction predicts the pressure dis- 
tribution quite closely. There is some evidence of 
deviations in the vicinity of the shoulder, but the pres- 
sure returns nearly to the free-stream static value as 
close as a half a radius aft of this point. The pressure 


distribution calculated from measurement of the 
Mach lines on the schlieren photograph (Fig. 5) is 


shown in Fig. 6. 


With these pressure distributions in mind, the heat- 
transfer measurements were confined to an area one 
or more diameters back from the shoulder in order to 
Both thin and 
calorimeter heat-transfer gages, as described in references 
The value of the transition 


attain a zero pressure gradient flow. 


11 and 15, have been used. 
Reynolds Number obtained by the two types of gages 


did not differ appreciably. 


The results of these measurements are shown in 
Fig. 7. Natural transition can be seen to occur at a 
Reynolds Number of approximately 4.0 to 6.0 * 10° 
based on the fluid properties at the edge of the boundary 
layer. Transition was obtained at a Reynolds Number 
of about 100,000 in the case where a circumferential 
trip wire was used near the sonic point. The Reynolds 
Number based on the momentum thickness was of the 
order 100 at the trip wire when transition occurred. 


It can be seen that the scatter of the data increases 
markedly at the high end of the Reynolds Number 
range. It has been found that, at these conditions, 
the combination of large heat-transfer rates and large 
viscous forces, as well as the extremely high impact 


+ An approximate analytical treatment for pressure gradient 
and body radius effects can be found in reference 16. 
t The diameter of the trip wire was of the order of the 


boundary-layer thickness. 


---= 


T TRANSFER ray, 


pressures, create a difficult problem of maintaining 
the heat-transfer gage integrity for a significant portion 
of the test time. As a result, there is a somewhat 
arbitrary selection of this data and it is conceivable 
that some of the scatter represents data which have 
not been fully corrected or interpreted for possible gage 
deformation. 

It can be noted in Fig. 7 that the data with the 
tripped boundary layers and with natural transition 
overlap in the Reynolds Number range covered. 
There is no evidence that the calorimeter gagest7 
impede the boundary-layer flow in a way to cause 
earlier transition or artificial thickening of the bound 
ary layer. 

In the transition range, using both tripped and 
untripped gages, it was found very difficult to analyze 
the records of the heat-transfer gage output. There 
is considerable evidence that the rates alternate from 
laniinar to turbulent values, but neither rate is clearly 


established. 


(b) Reynolds Number Dependence of Turbulent 

Measurements 

In addition to the scatter resulting from high Reyn 
olds Numbers, there is also some apparent scatter 
of the data caused by variation of the enthalpy ratio 
across the boundary layer associated with the different 
shock Mach Numbers used in the experiments. To 
minimize this effect, the data have been shown in 
Fig. S separated into experiments performed in three 
shock Mach Number ranges; 6-S at the high Reynolds 
Numbers, S—10, and 10-14 at the lower range of the 


tt The calorimeter heat-transfer gage consists of a strip of 
platinum approximately 0.001-in. thick, mounted directly on the 
surface of the model. Consequently, it protrudes into the 
boundary layer. The details of the calorimeter heat-transfer 


gage are discussed in references 11 and 15 
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Reynolds Number scale.j This graduation of the 
shock Mach Numbers at which experiments are per- 
formed at various Reynolds Numbers is unfortunately 
one of the operating characteristics of the shock tube. 
It is primarily the result of the structural limitations 
placed on the driver for safety reasons and also is due 
to attenuation in shock strength, particularly at the 
high Reynolds Numbers. 

It can be seen in Fig. 8 that, in all the cases tested, 
the value of the turbulent heat-transfer rate is reason- 
ably predicted by the incompressible turbulent skin-fric- 
tion coefficient and Reynolds analogy modified for 
Prandtl Number dependence, as given by Eqs. (33) and 
(34), with the Lewis Number set equal to unity. This 
trend shows that the incompressible skin-friction co- 
efficient is applicable to the highly cooled turbulent 
boundary layer, and, in particular, the Reynolds 
Number dependence (/e,)°-’ is verified. 

At the low Reynolds Number end, some of the high 
shock Mach Number data (see Fig. 8) tends to depart 
slightly from the eight-tenth power Reynolds Number 
dependence. This trend can probably be accounted 
for by the fact that the turbulent boundary layer is 
not completely developed at this low Reynolds Number. 

There appears to be some trend, indicated by the 
three curves in Fig. 8, for the measured heat-transfer 
rates to be greater than the incompressible result as 
the enthalpy ratio h,/h,, is increased—.e., higher shock 
Mach Numbers. 


(c) Effects Due to High-Temperature Fluid Properties 


In order to investigate the effect of variations of 
fluid properties, including the degree of dissociation, 
experiments were performed at several Reynolds 
Numbers over the range of free-stream-to-wall enthalpy 
ratios accessible. One set of experiments was _per- 
formed at a Reynolds Number, Re,, of 2.85 X 10° 
(plus or minus 10 per cent) using a wire at the sonic 
point for producing the turbulent boundary layer. 
Another set of experiments was carried out at a Reyn- 
olds Number Ke,, of 8.0 X 10°, for which natural 
transition existed. The free-stream-to-wall enthalpy 
ratios varied from 8.0 to 45.0 (shock Mach Number 
6.0 to 13.5). The highest degree of dissociation at- 
tained was a state in which 27 per cent of the original 
molecules were dissociated. 

The results from these experiments are shown in 
Fig.9. The data are normalized with respect to Reyn- 
olds Number in order to make correlations on the 
other parameters more easily discernible. The values 
of the normalized heat-transfer rates, g/Re,°-’, cal- 
culated for typical shock tube conditions from the 
incompressible skin-friction coefficient is shown by the 
solid line on this Figure. This line corresponds to 

t In order to establish a frame of reference, the variation of 
the other dependent flow parameters are as follows: 


WU he/hy Te/Tw M, 

5-8 8-16 7.0-11.5 1.70-1.95 

8-10 16-25 11.5-14.5 1.90—-2.10 
2.00—2.20 


10-14 25-50 14.5-23.0 


Eqs. (33) and (35) with the Lewis Number correction 
factor set equal to unity. 

The ordinate of an experimental point divided by 
the value of the ordinate of the curve at the same 
enthalpy ratio is just the ratio of actual to incompres- 
sible skin-friction coefficient, c;/c;;. It can be seen that 
the value of this ratio has a definite increasing trend 
with increasing enthalpy. The mean value varies from 
possibly slightly below the incompressible result at 
low cooling rates to approximately 1.20 at the highest 
enthalpy ratios. The results appear to scatter about 
a line through the mean of the data by approximately 
plus or minus 15 per cent. 

Methods of accounting for variations of c;/c,, with 
flow conditions have been suggested for low-tempera- 
ture flows. The reference temperature method of 
Rubesin and Johnson’? was interpreted by Eckert!® 
in terms of a reference enthalpy allowing extension to 
the conditions of these experiments. These methods 
determine the parameter 


92) 


(u* ue)?” (36) 


Cr Cy (p*/p.)*° 
where the starred quantities denote reference conditions 
which will reduce the compressible skin-friction coef- 
ficient to its incompressible form. 

The reference enthalpy—.e., the point in the bound- 
ary at which the starred properties are evaluated—is 


given by 


h* = h, + (1/2)(h, — h.) + 


0.22(0*)"3(H, — h.) (37) 


From the above expression, h*/h, has a value of ap- 
proximately 0.60 for the present experiments. 

The result obtained using the Eckert reference 
enthalpy method is shown in Fig. 9 as a dashed line. 
It can be seen that the prediction lies approximately 
20 per cent above the mean of the measured values of 
q/Re,°’. The increasing values of cy, c;, with larger 
enthalpy ratios obtained from Eckert’s method are 
primarily due to the effect of dissociation on the density 
ratio, p*‘p,. 

This trend of increasing values of c, ¢c,, with enthalpy 
ratio (degree of dissociation) can be correlated semi- 
empirically by accounting for the effect of the laminar 
Lewis Number, as discussed in Section (2c). In fact, 
if the form of the Lewis Number correction is taken to 
be given by Eq. (32), it is found that a value of 8 = 
1.0 gives a reasonable correlation of the data if the 
laminar Lewis Number is taken to be 1.4. It should 
be noted that the laminar Lewis Number is not ac- 
curately known. As a result, this correlation is not 
able to separate the effect of Lewis Number and @ in 
Eq. (32), but rather determines that the quantity 
(L* — 1) is approximately 0.40, 

In terms of Nusselt Number, the above correlation 
is represented by the following equation: 


Nu = 0.029 o/?Re,°3[1 + 0.4(hp,/H,) | (38) 


or, in terms of the ratio of the compressible to the in- 
compressible skin-friction coefficient, 
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ivi Jr = 1 + OA(hy, H,) (59) 


Cy/ Cy = 


In the above equation, the fraction of energy in dis- 
sociation can be calculated from the expression 


he = ho +hp, =~ €,T + hp, (40) 


For the case considered, the fraction of molecules 
dissociated is not greater than 30 per cent and the 
specific heat €, is very close to the value for molecules 
with full vibrational excitation 
where Jt) is the gas constant for air molecules. 

It was shown in reference 10 that, at free-stream 
Mach Numbers of the order of two, the value of 
cy Cy (for a free-stream-to-wall enthalpy ratio of 
unity) obtained from correlation of experimental data 
is approximately 0.95. The data of Fig. 9 does not 
appear to deviate from the incompressible result up 
to a free-stream-to-wall enthalpy ratio of 15, at which 
point dissociation first appears. Thus, it appears 
likely that the departure from the incompressible 
results is not due to changes in free-stream-to-wall 
enthalpy ratio, as this ratio has already changed by 
a large factor without significant effect on c,/c;;, but 
rather due to the dissociation of the air molecules as 
indicated by the proposed semiempirical correlation. 


(4) CONCLUSIONS 


The heat-transfer rate problem has been considered 
for the case of a highly cooled and partially dissociated 
turbulent boundary layer with the following results: 

(1) For the case of Prandtl and Lewis Numbers both 
unity, the equations of mean motion are essentially 
the same as for a nondissociating gas if the total 
enthalpy of the gas is considered to properly include 
the energy of dissociation. 

(2) From theoretical and physical arguments, if the 
Lewis and Prandtl Numbers are both equal to one, 
the incompressible heat-transfer result 
external fluid properties is expected to apply for the 
highly cooled turbulent boundary layer with zero 


based on 


pressure gradient. 

(3) Semiempirical correction factors are presented to 
account for values of the Prandtl and Lewis Number 
different from unity. 

Shock tube experiments have been performed and 
the results are discussed. It has been found that: 

(4) Transition occurs in a shock tube ona blunt body, 
and turbulent heat-transfer rates can be measured. 

(5) Turbulent heat-transfer measurements in the 
region of zero pressure gradient yield results close to 
predictions from the incompressible theory. 

(6) The deviations from the incompressible theory 
were found to be less than 15 per cent for the conditions 


SPACE 


ie., Cp ~ (9/2) Rol 
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covered by the experiments. The departure from the 
incompressible theory can be correlated by a correction 
factor, depending on the degree of dissociation of the 


gas and the laminar Lewis Number. 
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Theoretical Investigation of the Flow Field 
About Blunt-Nosed Bodies in Supersonic 


Flight’ 


ROBERTO VAGLIO-LAURIN*® ann ANTONIO FERRI** 
Polytechnic Institute of Brooklyn and Gruen Applied Science Laboratories, Inc. 


SUMMARY 


A numerical method for obtaining the solution to the inverse 
problem of the flow behind a given detached shock to any desired 
is presented. The cases of zero and small incidence 
The combination of sets of such solutions satis- 


accuracy 
are considered 
fying prescribed boundary conditions (body shapes) is described. 
Particular attention is devoted to the analysis of the sonic and 
subsonic region of the flow field. Convergence and stability of 
the stepwise integration from the shock in the elliptic region are 
discussed. Numerical examples are also included 


(1) INTRODUCTION 


AS A CONSEQUENCE of the large aerodynamic heat 
input in hypersonic flight, blunt shapes have be- 
come of interest for bodies and wings in contrast to the 
practice at conventional supersonic speeds. Thus, the 
aerodynamicist is faced again with problems of the 
transonie type, wherein regions of subsonic and super- 
sonic flow are encountered simultaneously. Consider- 
able effort has been devoted, by many authors, to the 


these problems. It that 


understanding of appears 
one can distinguish between two categories of flow 
fields: (a) flow fields where the domain of influence of 
the sonic line in the hyperbolic region does not include 
portions of the body surface [Fig. 1(a)]; and (b) flow 
fields where portions of the body surface fall within the 
afore-mentioned domain of influence [Fig. 1(b) ]. 

To a first approximation, the problems of category 
(a) can be analyzed by methods based on the assump- 
tion that the existence of a sonic line is of secondary 
importance (reference 1), while the problems of cate- 
gory (b) show an important influence of the shape of 
the body in the sonic region (reference 2). An analo- 
gous distinction applies when more accurate methods of 
analysis (numerical methods) are developed. Indeed, 
the sonic line in the flow field of Fig. 1(a) can be deter- 
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mined independently of the subsonic flow in the region 
ABDE when the portion AC of the shock is prescribed, 
while the determination of the sonic line in the flow 
about the body of Fig. 1(b) requires a simultaneous 
solution for the elliptic region ABDE. Additional 
complications arise when singularities are present in the 
flow field, as, for example, at the sharp corner D in 
Fig. 1(b); in such a case a preliminary analysis of the 
singularity is required before one may proceed to the 
matching of the subsonic and supersonic portions of the 
flow. 

The present paper describes a numerical method of 
analysis for the subsonic and transonic portions of the 
flow fields in question. The solution for the type of 
body shown in Fig. l(a) (‘rounded body’’) is obtained 
by an inverse procedure whereby a family of shock 
shapes AC are prescribed and the boundary conditions 
at the surface of the given body are satisfied by assum- 
ing a continuous variation of the flow properties with 
the shape of the shock, without linearization of the 
equations of motion. 

The solution for the type of body shown in Fig. 1(b) 
(“flat body’’) involves different sets of initial data and 
proceeds along different lines, according to whether the 
body profile is continuous or discontinuous at the 
shoulder. In the former case, considerations on the 
analytic continuation of the solution in the complex 
space for several representative shock shapes indicate 
that the flow field between the shock and the body is, 
in general, analytic. Under these circumstances, it is 
still convenient to seek the solution by an inverse 
method analogous to that used in the case of the 
different initial 


However, sets of 


namely, consistent shapes of 


“rounded body.” 
values are prescribed 
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the shock AB and of the sonic line BD.* With these 
data the solution in the elliptic region proceeds along 
the same lines, as in the previous case, and the flow 
about a given body is obtained by assuming a continu- 
ous variation of the flow properties with the prescribed 
boundary conditions. 

Among the many interesting cases of discontinuities 
in the body profile, attention has been devoted mainly 
to those characterized by a sharp corner which deter- 
mines the location of the sonic point. The correspond- 
ing sonic singularity can be investigated by extending 
the results established by Guderley*® for two-dimen- 
sional irrotational flow. If either a shape of the sonic 
line or a pressure distribution on the body are prescribed 
consistent with the afore-mentioned singularity, then 
the determination of the supersonic portion of the flow 
field included between the sonic line and the limiting 
characteristic can be performed by numerical methods. 
When the shape of the sonic line is given, the deter- 
mination of the elliptic flow field may proceed either 
from the shock (to be selected consistently) or from 
the body (where a consistent pressure distribution is 
prescribed), the final solution being obtained by iter- 
ation. However, if one intends to proceed from the 
body toward the shock, it is more suitable to follow 
the second alternative mentioned above —namely, one 
prescribes a pressure distribution consistent with the 
singularity at the shoulder and then analyzes simul- 
taneously the elliptic and the transonic regions up to 
the limiting characteristic, the final solution again being 
obtained by iteration. 

In all the cases considered, the method of analysis 
involves characteristic procedures for the supersonic 
portion of the flow field and a forward numerical inte- 
gration of the equations of motion in the elliptic region. 
The reliability of the latter technique has been investi- 
gated, both from the point of view of stability and of 
convergence of the solution; it appears that it can be 
quite satisfactory, and the available numerical results 
corroborate this conclusion. 

The various problems outlined above are treated 
here in detail along the following lines. The method 
of analysis for the flow field behind a given axisym- 
metrical detached shock is described first with par- 
ticular emphasis on the sonic and subsonic region. 
Convergence and stability of the solution in the elliptic 
domain are discussed. The combination of sets of 
such solutions satisfying prescribed boundary condi- 





* Observe that a ‘‘smoothing”’ process is included implicitly in 
this approach; strictly speaking, the shock and the sonic line can 
not be prescribed independently. Also observe that, for a given 
shape of the sonic line, no flow property there can be imposed 
arbitrarily; indeed, the momentum equation in the direction 
tangent to the sonic line, Bernoulli’s equation and Crocco’s 
theorem must be used in the determination of the velocity direc- 
tion, of the pressure and of the density. Corresponding to each 
set of boundary values one must determine the portion of body 
profile in the supersonic region which affects the sonic line. This 
portion of the profile must also be considered in the final combi- 
nation of the known solutions satisfying the prescribed body 


shape. 
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tions (body shapes) is then indicated. The flow about 
a body at small incidence is considered subsequently 
within the approximation that the three-dimensional 
effects may be linearized; the validity of this ap. 
proximation is discussed. Finally, the results of one 
complete numerical calculation are presented. 


(2) AXISYMMETRIC BoDy AT ZERO INCIDENCE 


The case of a ‘rounded body”’ will be discussed first 
[Fig. 1(a)]. As mentioned in the Introduction, the 
flow field about such a body is obtained by combining 
a suitable set of solutions, each corresponding to a pre- 
scribed axisymmetrical detached shock that com- 
pletely defines the sonic line and the elliptic region. 
Only the flow of a perfect gas will be considered here; 
however, the method can easily be extended to include 
real gas effects to any desired accuracy.‘ 

For a given shock shape one begins with the numeri- 
cal determination of the flow field in the supersonic re- 
gion BCD, applying the method of characteristics. 
A standard iteration procedure may be used, except in 
the neighborhood of the sonic line, where a particular 
modification must be introduced.’ This portion of the 
analysis is rather straightforward and does not require 
additional consideration. However, it may be of in- 
terest to notice that if the portion AC of the shock is 
described by one analytic law, no singularities should 
be encountered in the whole domain ACDE? and, thus, 
the integration of the equations governing the motion 
may proceed along the same network used in the sub- 
sonic region. Both procedures have been used in the 
numerical example discussed in Section (6); identical 
results were obtained. 

The distributions of flow properties immediately be- 
hind the shock AB and along the sonic line BD consti- 
tute the initial values for the forward integration of the 
equations governing the motion in the elliptic region 
ABDE. Setting aside the question of analytic con- 
tinuation of the solution in the complex domain, and 
the questionof stability (to be discussed in the following 
section), there appears no a priori reason why the calcu- 
lations should not be carried out in the physical plane. 
Practical experience, however, has shown that, if one 
tries to establish a Taylor series expansion for the en- 
tropy in the neighborhood of the sonic point on the 
shock, this expansion exhibits a very poor convergence 
at least up to the third term for some shock shapes— 
e.g., a parabolic shock. Thus, satisfactory accuracy 
can be obtained only at the expense of a large amount 
of numerical work. Therefore, it was found convenient 
to perform a change of independent variables from the 
physical plane x, 7 (where x coincides with the axis of 

7 Singularities in the supersonic region would appear at the 
shock; also, they would be reflected entirely at the sonic line in 
the same way as in plane potential flow, because local effects are 
two-dimensional and the entropy terms in the equations of char- 
acteristics vanish at the sonic velocity. Nonexistence of singu- 
larities in the elliptic region, when the shock is given by one 
analytic law, is discussed in Section (3) of this paper. 
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svmmetry of the body) to the plane 7, y where 


7T=2y/r*>;> y=r (1) 
and y is the stream function. In the plane 7, y the 
body is represented by the line 7 = 0, the shock by a 


line 7 constant, the axis of symmetry by the line 
y = 0, and the streamlines of the flow by the curves 
ry? = constant (Fig. 2). 

It is of interest that the transformation in question 
is useful also in constructing a procedure for the 
direct solution of the problem, based either on suc- 
cessive forward integrations beginning from assumed 
pressure distributions on the body, or on a relaxation 
procedure for the entire elliptic region. Indeed, a 
priori knowledge of the location of the streamlines and 
of the shock facilitates the guess of a consistent entropy 
distribution and the application of the boundary condi- 
tions at the shock front. 

In detail, the analysis proceeds as follows. 
the free-stream velocity, pressure, and 


Denote 
by V., Po» Pp 
density. In the region of disturbed flow indicate by 
V., U and V., V the x- and r-components of velocity re- 
ferred to a cylindrical coordinate system having the 
origin at the vertex of the shock (Fig. 1); also denote 
the pressure by p.. |’..°P, and the density by p.. X6 with 


6 = (vy — 1) /(¥ +1) (2) 


In this notation the equations of conservation of mass, 
momentum, and entropy are (for the axisymmetric 


case ) 
(0/Ox) (R Ur) + (0/0r) (R Vr) = 0 
U(oOU /ox) + V(OU /or) = —(6/R) (OP ox) | 3) 
U(OV dx) + V(OV/dr) = —(6/R) (OP °r)| ‘ 
U(0 0x) (PR~”) + V(0/0r) (PR~7) = 0 


The first equation in the system (3) implies the existence 
of a stream function W such that 


Oy /Ox = —(RV1r/5), Ov/Or = RUr/b (4) 


Upon introduction of the independent variables 7, y, 
defined by the Eqs. (1), the system of Eq. (3) can be 
transformed into the following: 
(0/07) [((U/V) — (76/RV)] + 
(6/2) (0/dy) (v/RV) = 0 
(0/dr) (P + rU) — 
[(1/2y) (0/dy)] (Uy?) = 0 (5) 
U? + V? + (1 + 4) (P/R) = 
1+ [((1 — 6)/6M,?] 
PR™* = f(ry?) 
where f(7y*) is an arbitrary function determined by the 
shape of the shock which maps onto the line 7 = 1. 
The first two equations in the system (5) are readily 
obtained from the first two equations in the system (3) 
by using the relations 
0/Oox = —(2RV/yés) (0/07) { : 
0/dr = (2/y) [(RU/5) — 7] (0/07) + (0/Oy) (6) 
The third equation in the system (5) is Bernoulli’s equa- 
tion; the fourth equation is a direct consequence of the 
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plane and its representation in the physical plane. 


conservation of entropy along the streamlines. The 
pressure and the density can be eliminated from the 


Eqs. (5) to obtain 


P [(c — U2 — V2)/(1 +6) fJo rr ft 
= = o\ ¢1(1 —8)/28 (7a) 
R [(e — U? — V*)/(1 + 4)f] f ; 
with c=1+ [((1 — 6)/1M,°6] (7b) 


and two simultaneous nonlinear partial differential 
equations for the two unknown velocity components 


\ 


A,(OU/Or) + B,(OV/Or) + 
C\(OU/dyv) + Di(OV/oy) + EF, . (8) 
A,(OU/dr) + BOV/dr) + 
C,(0OU/dy) + Ey = 0 


where 
A, = V}1 4+ [(6 — 1)/(6 + 1)Jrufo~* we 3 
[(¢ — 77 — TR/O 4. §)]~ fre 26)} (9a) 


A, = r — (U/S8) [(e — U2? —V?)/(1 + 8)f]" 7?” -— (9b) 


By = —U+7r[(c —U? —V?)/(1 + 5)f]o7 YY’ & 
{6 + [(6 — 1)/(6 + 1)]V? xX 
l(c — U2 — vy?)/(1 + 6)] . (9c) 


B, = —(V/8) [(e — U2 — V2)/(1 + df" 7?” = (9d) 
C, = —[(6 — 1)/2(1 + 4)]-(v/f) X 

UV[(c — U? — V2)/(1 + 6)f]7 OF?’ = (Ge) 

C. = —(y/2) (9f) 
D, = —(y/2) [(e —U? —V?)/(1 + &)f]°~? x 
16 + ((6 — 1)/(6 + 1)]V2[(e — U2? — V2) + 

(1 + 6)]-"} (9g) 

E, = —(V6/2) [(c — U? —V?)/(1 + 6)f]°~ = (Gh) 


Ey = [(6 — 1)/26] [((c — U2 — V*) + 
(1 + &)f]?r?’*. (of/Or) (9) 


At the axis (y = 0) Eqs. (8) become, respectively, 


oV/or = 0 ) 

Ul(ce — UX/(1 + afo]"-2/ = wf 
Thus, U at the axis can be computed without solving 
for the entire flow field. However, the inverse trans- 
formation to the physical plane (x, 7), which, in general, 
entails the integration of the equation 


dx = —(1/RV) [(y6/2)dr + (76 — RU) dy] (11) 


and, in particular, the integration of the equation 
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dx = —(U/2r) (OV/dy)—'! dr 


at the axis, can only be performed upon completion of 
the flow field calculation. 

The boundary conditions at the shock front described 
by the equation x = F(y) are 


U=1-(1-—6)) [1d + F)) —-(d W.°)} 
V= (1-6) F, [1d + F,)] - (/1.°)! 
f(ry?) = (1 — 6) ) fl + F,2)] - 

(6/(1 + 6) M.2]}} (1 + [UA — 6) 617.2] X 

(1 + F,2}0+H/a-8 

where F, stands for (d‘dy)F(y). Since the shock 
transforms into the line 7 = 1 in the plane (7, y), the 
Eqs. (12) also permit the determination of the partial 
derivatives (O” Oy”) of the flow properties at the shock 
itself. For example, 


oU/dy = 2(1 — 6) [F, F,,/(1 + F,?)?] j 
OV /oy = (1 — 5)F,, X 

fa — F2)/ + £2] - (1/2) 
Upon evaluation of Eqs. (12) and (13), Eqs. (S) can be 
solved algebraically for the derivatives Ol” O07 and 
OV Or. 


by a simple extension of this procedure. 


Higher order derivatives could be determined 


The shape of the sonic line in the 7, y plane can be 
obtained by a point to point transformation from the 
physical plane where the velocity components and the 


entropy are known. Indeed, it is 


1 — 6) AS 


or [(1 + 6)/(1 — 6)]6"*” Mn? f(ry?) =e 


and y = r. Upon solution of this equation for f(ry?), 
7 at each point can be determined from the third of 
Eqs. (12). 

Finally, the distribution of U’ at the axis can be ob- 
tained by solving Eq. (10). 

Given a shock shape, the flow field in the elliptic re- 
gion is determined by numerical methods, the system 
of differential Eqs. (8) being replaced by a system of 
difference equations. The solution proceeds step by 
step from the shock (line 7 = 1), each step resulting in 
the simultaneous evaluation of the flow properties at 
a fixed number of points on successive lines 7 = con- 
stant (Fig. 2). 
integration are required to satisfy the prescribed values 


In each step, the results of the forward 
on the sonic line and on the axis. The first step, which 
results in the determination of the velocity components 
U and V at a certain number points on the line O—p 
(Fig. 2), is described in detail here; the same procedure 
can be repeated successively from the afore-mentioned 
line on to the body (line r = 0). 

Indicate by the general point on the line O—p and 
by us the point having coordinates (1, y,). Eqs. (8) 
applied at this general point can be expanded into 
finite difference form, retaining second-order terms, by 


substituting 
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(OU /Or), = 2[(Un — Uns)/(ta — tns)] - 
(OU Or), (l4a 
‘ica Uy Ota  ¥ i) | 


(OU/Oy), = [(Uns 


for the derivatives, and by approximating U’ and J’ in 
the expressions of the coefficients A; through /» [Eqs 
(9a) through (91) ] with expressions like 


U, = U, 


+ (OU'/Or)n: (Tr — Tns) (14b) 


After the initial values at the shock front have been 
computed from Eqs. (12) and (13), Eq. (14b) permits 
the evaluation of the coefficients A; through /. to the 
desired approximation; subsequent substitution of 
Eqs. (14a) into Eqs. (8) give place to two linear equa- 


tions that can be put in the following form: 


—CinUn-1 + a:,U, + CnrU ens — dinV, iT 


DinVn t+ dinVnar + Ci, Q (15) 
where 7 = 1,2;1<" < p — lj; and 
lin = BA, b,, = 2B, 
Cin = Cin[(ta — Tns)/(Yn41 — Vn-1)] 
din = Din{(tn — Tns)/(Ynt1 — Vn—1)] P 
_ > ( ) 
Sn = ~HA NM + BaVnl — ” 
[A ;,(CU' O7),, + Bi,(OV Or),. — 
Ein] (tn — Tne) 


with the convention D,, = 0. The unknowns Ul’, and 


Vv, 1 Sa" < p — 1) are determined by solving the 
system of 2(p — 1) linear equations obtained from Eq. 
(15) by setting 7 = 1, 2 and = 1, 2, (p — 1). 


These values of U’,, and J”, are then used to evaluate the 
velocity derivatives (0 Oy), and (O O7r),; thus, one can 
proceed to the determination of the flow properties on 
the next line 7 = constant, and so on. 

A few additional comments are in order concerning 
the calculation of the flow field about a ‘‘flat body” 
(Fig. Ib). As stated in the Introduction, analytic 
shapes for the shock AB and the sonic line BD are 
selected in connection with a prescribed body having 
an analytic profile. Initial values immediately behind 
the shock are determined from Eqs. (12) and (13) pre- 
sented above; initial values along the given sonic line 
are obtained consistent with its shape* by numerical 
integration of the momentum equation in the direction 
tangent to the sonic line combined with Bernoulli's 
equation and with the equation expressing conservation 
of entropy along the streamlines. The solution in the 
elliptic region then proceeds along the same lines as 
indicated above. The portion of body profile in the 
supersonic region which affects the sonic line may be 
determined by characteristic procedures.® 

In the case where the profile of the body exhibits a 
sharp corner and thus determines the location of the 
sonic point, the pertinent singularity can be obtained 
from that occurring in plane potential flow by adding 
analytical correction terms. This conclusion can easily 
be reached if one observes that (1) a centered expansion 


* The shape of the sonic line must be consistent with the shape 
of the shock to the same order of accuracy (order of the deriva- 
tives) carried in the numerical calculation. 
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alwavs exhibits a two-dimensional behavior in the im- 
mediate neighborhood of the corner, and (2) the changes 
of flow properties on the singular streamline (body sur- 
face) are governed by the equation valid along a char- 


acteristic line—namely, ! 


cot u(dg g) -- dé + sin w cos uw (dS yR) = 0 (17) 


and, thus, are not affected by the entropy gradient in 
the direction normal to the body surface. Hence, the 
shape of the sonic line and the pressure distribution on 
the body in the neighborhood of the corner can be 
analyzed following a method analogous to that described 
in reference 6 (Guderley’s results apply directly only 
in the case of small perturbations). With these data, 
the solution can be continued by numerical methods 
following the lines described in the Introduction. Ex- 
cept for the modifications just described, which are 
connected with the determination of the initial values, 
no changes are required for the numerical determination 
of the subsonic portion of the flow field. 

The question of stability, which is frequently raised 
for the procedure proposed here, will be considered in 


the following section. 


(3) STABILITY AND CONVERGENCE OF THE SOLUTION 
IN THE SUBSONIC PORTION OF THE FLOW FIELD 


Hadamard’ first pointed out examples of instability 
for solutions of the Laplace equation satisfying given 
initial values on an open boundary; in particular, he 
showed that the solution may depend on the initial 
values in a discontinuous fashion. For this reason, the 
problem in question has been called an “incorrectly 
posed Cauchy's problem.’’ However, there has been 
much speculation recently as to the real bearing of this 
conclusion. An interesting result in this connection 
has been obtained by Pucci,’ who has proved the fol- 
lowing theorem. Consider a function L(x, y) harmonic 
in a simply connected domain, which includes the in- 
terval (0, 1) of the x axis (Fig. 3). Let / be the dis- 
tance between the segment (0, 1) and the boundary 
§D. For arbitrary values of the index », give approxi- 
mate values a,” and 6; for the function U(x, y) and 
for its derivative U’,(x, y) at the (42 + 1) points (7/47, 
0),7 = 0, 1,...4n, of the interval in question, such 
that 


U(i/4n, 0) — a; | < E,; 


U,(t/4n, 0) — BY” | < EB, 


Corresponding to this approximation define the func- 
tion U’,(x, vy) 


U(x, vy) = do (-1)* [y?* (2s) JA a yy) + 


s=0 


> (—1)*[y*+!/ (2s + 1)!JA* Bi) (18) 


0 
where the operator A stands for 
Aa,” = 4n(a™ ini — a;”) wheni < 2n 


Aa;™ = 4n (a; — a1) wheni> 2n 
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Fic. 3. Domain considered in Pucci’s theorem 


If the error £,, of the approximate initial values satisfies 
the condition 

lim (£,,e°"") 0 (19) 

the succession ; U’,(x, y){ converges uniformly to U(x, 

y) in the rectangle R=(OSx<1,/y 

h’ is any positive number smaller than /. 

The importance of the theorem lies in the fact that it 


< h’) where 


indicates the possibility of constructing convergent 
numerical approximations to the solution according to 
Eq. (18); in the limit of m —~ © this coincides with the 
Taylor expansion of the function L(x, y). 

From the result in question one may conclude that 
satisfactory solutions of initial values problems can be 
obtained also for elliptic equations more general than 
the Laplace equation if a suitable pattern of integration 
The integration can proceed step by step 


0,0 <x < 1) to 


is followed. 
from the initial line (say, line y 
successive lines (y = constant, 0 < x < 1), provided 
that the flow properties at all the (42 + 1) points con- 
sidered in any particular step are determined within the 
approximation required by Eq. (19). It is of interest 
to notice in this connection that the degree of accuracy 
increases exponentially with the number of points 
carried* and with the magnitude of the step. With 
these two parameters fixed, the calculation must be set 
up in such a way as to maintain, in each step, the same 
approximation (errors beyond a certain decimal place). 
It appears that the truncation error can be controlled 
by using the difference approximation indicated in 
Pucci’s theorem, and that the number of terms to be 
carried in Eq. (18) is determined by the condition that 
the rounding error satisfies. Eq. (19) to the desired 
number of decimal places. If the condition on the 
rounding error requires that m terms be kept in this 
series, the numerical work can be facilitated by solving 
* This dependence is in qualitative agreement with the sug- 
who pointed out the possibility 


) 


gestion of Zlotnick and Newman,’ 
of obtaining satisfactory results by filtering high frequencies 
from the solution. 
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simultaneously for all the points involved in each step; 
indeed, with this approach only the first (7 — 1) terms 
must be evaluated on the initial line while the mth term 
results from the solution of the pertinent system of equa- 
tions. For example, the method of analysis presented 
in the previous section considers second-order terms in 
each step, even though only first derivatives are evalu- 
ated on each initial line. Since Eq. (18) is a finite 
difference Taylor expansion of the unknown function 
U(x, y), the modification in question does not alter 
the difference approximation within the selected order 
of accuracy (7 terms). 

Whatever the procedure followed in a practical cal- 
culation, one must ascertain the accuracy of the results 
obtained in each step. This check can be performed a 
posteriori in the following way. Upon completion of 
the calculations relative to a certain step where uth 
order derivatives have been considered, evaluate the 
derivatives of the (7 + 1) order by finite differences and 
assess whether their contribution to the results is 
negligible, consistent with the requirement of Eq. 
(19). If this condition is not satisfied, the calculation 
must be repeated, for example, by reducing the magni- 
tude of the step; only a small amount of additional 
numerical work is then involved. 

In the light of this discussion it appears that a study 
of the propagation of errors should also be carried out 
according to the provisions discussed above; indeed, 
if the requirement of Eq. (19) is violated, the accuracy 
of the solution is open to question. 

In view of the assessed possibility of constructing 
stable solutions within a rectangular region having the 
initial line as basis and one vertex on the boundary of 
the domain of analyticity of the solution, it is of interest 
to investigate the extension of the latter domain. The 
following question will be considered in this connection: 
given an analytic shape of the shock, what is the radius 
of convergence of a Taylor series expansion for the 
solution, initiating at the shock? The answer can 
easily be obtained by continuing the solution in the 
complex domain and by investigating the possible ap- 
pearance of singularities. Some results of Garabedian” 
will be used in the following considerations. 


Garabedian formulates the Cauchy problem for the 
determination of the flow field behind a given detached 
shock in the complex domain.* For this purpose he 
transforms the elliptic partial differential equation for 
the stream function in the subsonic region into a canoni- 
cal hyperbolic system of equations by taking combi- 
nations = and ¥ of the complex characteristic variables 


aand pg 
= (a+68)/2, ¥ = (a — B)/2i 


as independent variables and, subsequently, by keeping 
& real while letting ¥ = — + 7% be complex. Upon this 


* This is equivalent to introducing a time dependence for the 





solution of the problem. 
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classical sense.f With the independent variables 
Z, &, and 4 the shock can be mapped into a portion of 
the plane ¢ = 0 while the characteristics in a plane £ 


= constant are defined by the equations * + 4 = con- 


stant and  — 7 = constant. If the actual shock is 
described by an analytic law, the initial values in the 
plane ¢ = O are also analytic, and no singularity jis 


carried along the characteristics that intersect the initial 
curve in the complex domain. Singularities can only 
be carried by characteristics that either are tangent or 
do not intersect the aforesaid curve. For a given shape 
of the shock one can identify the point through which 
the tangent characteristics pass; however, the entire 
solution must be determined before one can perform 
the inverse transformation from the space of the vari- 
ables «, ¥ to the physical plane. A preliminary assess- 
ment is possible only if some simplifications are intro- 
duced; for example, one may assume that the equation 
under consideration has the same fixed complex char- 
acteristics as Laplace’s equation. 

Within this approximation, however, the desired re- 
sult can be obtained more readily by considering the 
complex space x, &, n, defined in terms of the actual co- 
ordinates x and r by the relations 


r=é&+17n 


x = X, 


With this change of variables the Laplace equation 
transforms into the wave equation (thus, the fixed 
slope of the characteristics in any plane £ = constant 
is known). The intersections between the surface 
obtained by analytic continuation of the law describing 
the shock and the planes € = constant can easily be 
determined. Then one need only to consider the char- 
acteristics tangent to the afore-mentioned intersections; 
these lines fix the range on the x axis where a Taylor 
series based on the given initial data will certainly con- 
verge. For example, consider a shock described by 
the hyperbola 


(x?/a") — (r?/b?) = 1 (20) 

The analytic continuation of this curve is 
(x?/a?) — [(E + i )?/b?] = 1 

and its intersection with the plane = 0 is the ellipse 

(x?/a?) + (n?/b?) = 1 (21) 
The characteristics of the wave equation tangent to the 
ellipse (21) are described by the equation 

"=x + Va? + b 


These points are the foci of the hyperbola described by 
Eq. (20). Hence, one obtains the result (also derived 
by Garabedian) that in the case of a shock having the 
shape of a hyperbola the first singularity on the axis 
would occur at the focus, if the actual equation govern- 
ing the flow were to have the same fixed complex char 


7 It must be noticed that with Garabedian’s method a very 
special treatment is required in the neighborhood of the sonic 
line where the equation of motion changes from elliptic to hyper 
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acteristics as Laplace’s equation. These considerations 
can be repeated for various planes § = constant; thus 
the loci of possible singular points (limiting lines) up- 
stream and downstream of the shock can be determined 
approximately and their location with respect to the 
body may be ascertained. 

In conclusion, we have shown that the numerical 
procedure proposed in this paper satisfies the require- 
ments to be inferred from Pucci’s theorem and that its 
range of convergence can be expected to cover the en- 
tire elliptic flow field at least for simple shapes of the 
shock. In any case, an approximation for this range 
of convergence can be determined a priori; if necessary, 
a more exact assessment can be performed while pro- 
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satisfied by perturbing the prescribed shape of the shock 
as a function of a parameter and by assuming contin- 
uous variations of the flow properties with the pertur- 
bation, without linearization of the equations of motion; 
(b) the boundary conditions are satisfied following 
a standard perturbation scheme which involves lineari- 
zation of the equations of motion. 

In either case the problem appears to be correctly 
set, and, thus, the solution depends continuously on 
the boundary conditions. The approach described 
at (b) is more rapid, particularly if the linearization is 
performed after transformation to the 7, y plane‘; 
however, the applicability of this approximation is 
limited to cases where the flow properties vary grad- 
ually throughout the flow field and particularly in the 








ceeding step by step with the calculation. On the basis 


of these considerations, satisfactory results can be ex- In other cases—e.g., 


neighborhood of the sonic line. 
flat bodies with shoulders having a large curvature 
the approach described at (a) is more rigorous. 

As for the numerical analysis, the aiternative (a) 
involves sets of standard axisymmetric solutions [cf. 
Section (3)], all proceeding in the same way from each 


pected from the proposed method of analysis; the nu- 
merical results presented in Section (6) corroborate this 
conclusion. 

The discussion of the zero incidence problem is con- 
cluded in the following section where the combination 
of sets of inverse problem solutions satisfying prescribed 
boundary conditions (body shapes) is considered. 


pertinent set of initial values; the alternative (b) in- 
volves a standard perturbation scheme. 

The discussion of the zero incidence problem is thus 
concluded. The flow about a blunt body at small inci- 


(4) DETERMINATION OF THE FLOW FIELD ABOUT A 
dence will be considered in the next section. 


GIVEN Bopy 


The determination of the flow field about a given 
body involves, in general, an iterative process or an 
equivalent perturbation scheme. Correspondingly, two 
alternative approaches may be followed: (a) the 
boundary conditions at the surface of the body are 


(5) FLow Asout A BLUNT Bopy AT SMALL INCIDENCE 


With reference to a cylindrical-coordinate system 
(x, r, 3), the equations governing the three-dimensional 
flow of an inviscid compressible gas are (where an 
asterisk indicates actual physical quantities) 


u*(Ou*/Ox) + v*(Ou*/Or) + w*(Ou*/rdd) = —(1/p*) (Op*/dx) 
u*(Ov*/Ox) + v*(Ov*/Or) + w*(dv*/rdd) — (w*?/r) = —(1/p*) (Op*/Or) 
u*(Ow*/Ox) + v*(Ow* /Or) + w*(Ow*/rdd) + (v*w*/r) = —(1/p*) (Op*/rod) (22) 


[O(p*u*)/Ox] + [0(p*v*r) /Or] + [0(p*w*)/rdd] = O 
u*(0/Ox) (p*p*~*) + v*(0/Or) (p*p* *) + w*(0/rdd) (p*p*~”*) = O 

For small incidence angles a let 
u*=V,(U+u) v®* =Vi.(V+), wt*=V.w, p* =p. 


where U’, V, P, R represent the solution for the flow at zero incidence and the quantities u, v, w, p, p represent small 
superposed disturbances due to the incidence; if one assumes that these disturbances can be linearized, the Eqs. (22) 


can be reduced to the following form: 


(0/dx) [Ru + pU] + (0/rdr) [rRv + rpV] + (R/r) (Ow/dd) = O 


u(O0U/dx) + U(du/dx) + v(OU/dr) + V(Ou/dr) = —(6/R) (Op/Ox) + (6p/R*) (OP, Ox) 
u(O0V/dx) + U(dv/dx) + v(0V/Or) + V(0v/Or) = —(6/R) (Op/Or) + (6p/R?) (OP /Or) (23) 
U(Ow/d0x) + V(Ow/dr) + V(w/r) = —(6/R) (Op/rdd) 


PR~* [U(Q/dx) + V(0/dr)] [(p/P) — y(e/R)] + [u(@/dx) + 2(0/dr)] [PR~?] = 0 


If the x axis is chosen coincident with the axis of symmetry, the boundary conditions at the surface of the body de- 
scribed by the equations x = X (r) are 


u/v = (dX /dr) (24) 
Within the approximation of linearization, the boundary conditions at the actual shock described by the equation 
x = F(r) + s(r, 8) can be reduced to conditions at the axisymmetric shock front; if the plane } = O is chosen coin- 


cident with the windward half of the pitch plane, the following relations are obtained on x = F(r): 
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u + (OU/dx)s = [4/(y + 1)] sin® B cos B(Os/Or) — [2/(y + 1)] sin B cos B [1 + (1/1/.? sin? B)] a cosd 
v+ (OV/Odx)s = —[2/(y + 1)] sin? B [cos? B — sin? B + (1/\/.2 sin? B)] (Os/dr) — 

[2/(y + 1)] cos? B [1 + (1/.A/..2 sin? B)]} @ cos J a 
w= asind + [2/(7y + 1)] sin? B[1 — (1/17. sin? B)] (Os/rdvd) -_ 
1/R? [p + (OR/Ox)s| = —[4/(y — 1)] (cos B/M..” sin® B) [sin? B(Os/07) — a cos J 
p + (0P/dx)s = —[4/(y + 1)] sin B cos B [sin? B (Os/Or) — a cos #] 


Inspection of Eq. 


turbances u, v, w, p, p, s can be factored out by taking u, 


, , ail 
u=au cos?, v= av costd, w=aw 


3) and (25 


Then Eqs. (2 


(0/dor) [ry Rv’ + r p’V] + (Rw’/r) 


u'(OU Ox) + U(ou'/ox) + peed Or) + V(ou’ 
u'(OV/Ox) + U(dv'/dx) + v’ (OV/dr) + V 
cata Ox) + V(Ow’/Or) + V(w’/r) = (6/R) 


PR~* [U(0/ox) + V(0/dr)] [(p’/P) — vy(p’, 


Ox)s’ = [4/(y + 1)] sin® B cos B(ds’/dr) — 


‘= —[2 


u’ + (OU 
v’ + (OV /Ox)s 


w’ = 1 — [2/(y + 1)] sin’ B [1 — (1/M..2 sin? B)] (s’ 
t’ + (OP/dx)s’ = —[4/(v + 1)] sin B cos B [sin? B(ds’ 
(1/R2) [p’ + (dR/dx)s'] = —[4/(y — 1)] (cos B/M..? 


It is of interest to notice that two particular solutions 
the system of Eqs. (27) are given by* 


=, 


oO 


u' = 0U/dr; v' = OV/Or; w’ = —(V/r) 


p’ = OP/dr; p’ = OR/Or; ds’/dr = (29) 
—(dB/dr) + [(y + 1)/2] [.2/(7.2 -— 1)] 
and by 
u’ = V+ x(0U/dr) — r(OU/Ox) 
y = ale x(OV/or) — r(0V/Ox) 
w= Vie). as jar = 1 (30) 
p’ = (OP ar) — r(OP/Odx), | 
p’ = x(OR/dr) — r(OR/Ox) 


Both solutions satisfy also the boundary conditions (24) 
and (28) in the neighborhood of the axis, when s’(0) = 
0 and the origin of the coordinate system is chosen coin- 
cident with the vertex of the axisymmetric shock. 
Since the solution (29) gives the perturbations pertain- 
ing to a small displacement of the body in the r direc- 
tion, while the solution (30) gives the perturbations 


* Professor G. B. 
attention of the authors. 


Whitham brought these solutions to the 
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Schematic diagram of the flow field and pressure 


distribution on the body. 


(23) and of the boundary conditions (24) 


sin d, 
) become, respectively (8 = cot 


+ (0/dx) [Ru’ 


(Ov’/Or) = 


(y + 1)] sin? 8 [cos? B — sin? B+ (1/M 


and (25) indicates that the dependence on # of the dis- 


p,S ~ acos J and w ~ a sin 3; namely, 


v, P, 


b= ap’ cosd, = ap’ cos? s=as’cosd’ (26) 
, 


| F,), 
+p Ul)= 


—(6/R) (Op’ 
—(6/R) (Op’ 


Ox) + 6(p’/R*) (OP /Ox) 
Or) + 6(p’/R?*) (OP /dr) (27) 


or) = 


(p’ r) 


‘R)] 4 [u’(0/Ox) + v'(0/dr)] [PR~7] = 0 


[2/(y + 1)] sin B cos B [1 + (1/M,,? sin? B) 
2 sin? B)] (ds’ 


11 — [2/(y + 1)] cos? 6 [1 + (1/M 


dr) — 


.2 sin? B)]}{ 


(28) 
r) 
dr) — 1] 


sin’ B) [sin? B(ds’/dr) — 1] 


pertaining to a small rotation of the coordinate system, 
the behavior of the actual solution in the neighborhood 
of the axis can be obtained by linear combination of 
(29) and (30); the determination of the coefficients in 
this combination is described below. 

Naturally, some discretion must be used in applying 
the linearized approximation, depending on the shape 
of the particular body under consideration. It is ex- 
pected that satisfactory results may be obtained for 
bodies having a quasi-spherical shape, while some 
doubt may be cast on the validity of this approach for 
flat bodies. It is also of interest to institute a compari- 
son between the subject boundary value problem and 
the problem encountered in the analysis of perturba- 
the transonic flow about a two-dimensional 
airfoil, The work of Cathleen Morawetz!! gives rea- 
sonable evidence that the perturbation problem for a 
blunt nose is correctly set; then the solution exhibits 


tions of 


a continuous dependence on the boundary conditions. 
To determine linearized perturbations belonging to a 


given body shape, one first computes “unit solutions” 
originating from prescribed simple perturbation of the 
shock shape (e.g., 8’n41 = Ar*"*!); subsequently these 
are superimposed to fulfill the boundary conditions ai 
the body. The behavior of each “ at the 
axis is obtained by combining the particular solutions 
(29) and (30) so as to match the prescribed ds’ /dr at 
the axis. After the initial values at the shock and the 
boundary values at the axis have been computed, the 


“unit solutions” is carried out by 


unit solution”’ 


determination of the 
forward integration in a fashion analogous to that de- 
tailed for the axisymmetric case. In particular, the in- 
tegration in the elliptic portion of the flow field is car- 
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ried out again in the plane of the transformed variables 7, y. With this transformation of variables Eqs. (27) be 
come, respectively, ; 


(1/y) (0/dy) py2[(2’/V) + (p’/R)]\ + (w’/V) + 2(0/d7) {Z — r[(e’/V) + (p’/R)]} = 0 
(0/07) (p’ + ru’) — (1/2y) (0/Oy) (u’y?) — [(v’/V) + (p'/R)] (OP Or) — Z(OU Or) = O 
Uu' + Vo’ + [C1 + 6)/2 [(p’ P) — (p'P/R?)| = 0 

(1/2y?) (0/Oy) (w’y*) — (0/07) (w’r) = (6/2RV)p’ 


a) p’ 1+ 6 p’ 1 Oo 6 1+ 6’ Z Of nuta/a—i 
- _ — y- — 7 R = 
Or ad 1—6R 2vovL° \P 1—6R P Or 


where Z = (R/5) [((U/V)v’ — u’ 32 





while the boundary conditions (24) and (28) remain formally unchanged except for a straight substitution of the 
variable y in place of the variable 7. The first, second, fourth, and fifth equation in the system (351) are obtained 
directly from the homologous equations in the system (27) by applying the relations (6); the 7 momentum equation 
[third in the system (27) | is replaced by Bernoulli’s equation [third in the system (31) |, which is a consequence of the 
system (27). Upon elimination of one unknown (e.g., p’) from Bernoulli’s equation, the remaining four equations are 
rewritten in the following form, which is convenient for the purpose of numerical computations :* 


A, (Ou’/Or) + As (Ov'/Or) + A2™ (Ow’/O7) + Ag (OP’/Or) + As“ (Ou’/OY) + 
A,” (Ov’/Oy) + Ax” (Ow’/Oy) + As\” (Op’/Oy) + B 0 (33) 


The coefficients A; and B“ are given by the following equations: 


Ay? = —2R{ (1/6) + [27U/(1 + 6)P]}, Ax = —(2/V) [r — (UR/8)], { 
A;Y = 0, Ay’? = —27(R/P?), (54) 
A; = [2y/(1 + 6)] (UR/P), Ac = (y/V), Ax? = 0, As? = y(R/P?) 
A,\® = 7, A, = A;% =0, As = 1, As® = —(y/2), Ac’ = Ary = As = 0 (35) 
A,® = A» = 0, A3© = —?2r, A, “= A; = Ag = 0, A; 3s As y= 0 (30) 
A, = —[2r/(1 — 6)] (UR/P), Ax® = A; = 0, A, = —(r/P) ate + 6)/(1 — 6)] (R/P) — lj, | 2 
A, = [y/(1 — 6)] (UR/P), A¢® = Az = 0, As = (y/2P) [1 + 6)/(1 — 6)] (R/P) — 145 : 
2 | 9) (=) O (\) eed ; E 7) (*) a) («.) 
BY = y — 2r — u’ + -}—2 -) + 
1+6L Ov\P or \ P 6 OF y ov \l Or \I 
>? (FF) —e | a) (5) a) (| )| ; 
2 > +TIy —2 b 
dr \eV VL oy \P2 "Or \P2 
Rou 26 UoP] , 1[UROU OP) , ROP , 
BO = - _ =i — + -— — p 
6 LOr 1+ 6POr ] 6 Or Or P? Or - (38) 


B® = w’ — (p’6/RV) 

a l k fe) ( ) fe) ( =) 1—d6R =| ; 1URoOf , 

Bo = y? —2 T + - — In" — : v + 
1—sdLyvyovy\ P Or a 6 f Or 6 Vf Or 


jl °F G8 .) ak ( +6R i) 
lQvy ov LP \1 —6P or LP \l1—6P i 


Eqs. (33) can then be applied at each point of the net- 


work used in the axisymmetric calculation, and there 
° me om . Sj ‘e >» enefficie : j > ev . P Eline 2 = 

be expanded in finite difference form by the same pro- Since the coe ficients of the syRem oO Eqs. 31) de 
cedure used for the determination of the basic flow. pend only on the basic flow, its characteristics a the 
There results for each line r = constant a svstem of complex domain will coincide with those of the basic 
. . . . : “ged _ > ( ‘ ; 7 ‘ < rti 4 y 7 > > ¢ 

4(p — 1) simultaneous linear equations in the 4(p — 1) flow; thus, the lomain of analy ticity of the pe rturba 
unknowns representing the values of the disturbances tion solutions is the same as that of the basic solution. 


, , 


u',v', w’, p’ at the (p — 1) nodal points on the line under Also, stability of the forward integration for the solu- 


consideration. Proceeding step by step from the shock tions in question can be obtained if the calculations are 
(r 1) to the body (7 = 0) the perturbation can thus carried out along the lines discussed in Section (3). 
be evaluated in the entire flow field. ; 
(6) NUMERICAL RESULTS 
“The superscripts (1) through (4) refer, respectively, to jai ; ‘ : ; 
coefficients in the continuity equation, in the x momentum Phe subject method - analysis has been applied = 
equation, in the 9 momentum equation, and in the energy equa- the determination of the flow field behind the axisym 


tion. metric shock x = r*/2. All the computations have 
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Fic. 5. Distribution of flow properties at various radial sections 
across the flow field. 


been carried out on a desk calculator; about 40 man- 
hours were required for the determination of the sonic 
line and about 100 man-hours were necessary to com- 
pute the whole elliptic flow field and to perform the 
inverse transformation from the plane (7, y) to the plane 
(x, 7). The interval 0 < 7 < 1 was divided in seven 
steps, and six to eight points* were carried in each step, 
as indicated in the transformed network of Fig. 2. 

The results, which are summarized in the plots of 
Fig. 4 and of Fig. 5, compare very satisfactorily with the 
data of reference 12. 


(7) CONCLUSIONS 


Numerical methods for obtaining tke flow field about 
blunt axisymmetrical bodies of arbitrary shape at zero 
incidence have been discussed. A method for the deter- 
mination of the effects of incidence within linearized 
approximation has also been presented. Attention 
has been devoted mainly to the analysis of the subsonic 
and sonic regions. The method, which is common to 
the axisymmetric and to the three-dimensional prob- 
lem, proceeds by forward integration from the pre- 
scribed values on the given initial curve, both in the 
hyperbolic and in the elliptic region. Criteria for in- 
suring stability of the solution in the latter portion of 
the flow field have been established, and the compliance 


* Points were added from the sonic line as the range in y in- 
creased while proceeding from the shock to the body. 
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of the proposed method with these criteria has been 
ascertained. On this basis it has been shown that the 
solution can be obtained rather rapidly even in the 
elliptic domain; a suitable change of variables that par- 
ticularly facilitates the analysis has been presented, and 
the pertinent procedure has been described in detail, 
A complete example of calculations has been quoted. 
The results are very satisfactory. 

The method in question also permits the investiga- 
tion of the flow field about bodies having singular pro- 
files—e.g., the ‘‘flat body”’ of Fig. 1. 

Beyond the interest connected with obtaining the 
complete transonic flow field about a blunt body, the 
results of the present investigation point to very prom- 
ising potentialities of suitable methods of forward 
integration for the determination of elliptic flow fields. 


REFERENCES 


1 Chester, W., Supersonic Flow Past a Bluff Body with a De- 
tached Shock—Part I—Two-Dimensional Body, J. Fluid Mech., 
Vol. I, Part 4, pp. 353-365, October, 1956. Part II—Axisymmet- 
vical Body, J. Fluid Mech., Vol. I, Part 5, pp. 490-496, November, 
1956. 

2 Hayes, W. D., Some Aspects of Hypersonic Flow, The Ramo- 
Wooldridge Co., Los Angeles, January, 1955. 

3’ Guderley, K. G., Singularities at the Sonic Velocity, WADC 
Rep. F-TR-1171-ND, June, 1948. 

4 Vaglio-Laurin, R., On the Analysis of the Flow-Field About 
Blunt-Nosed Bodies in Supersonic Flight, PIBAL Rep. 431 (to be 
published). 

5 Vaglio-Llaurin, R., Supersonic Flow-Field About Blunt-Nosed 
Bodies of Revolution—Part 3—The Methods of Analysis for the 
Subsonic and Sonic Regions, Gruen Appiied Science Labs., Tech 
Rep. No. 4, October, 1956. 

6 Vicenti, W. G., Wagoner, C. B., and Fisher, N. H., Calcula 
tions of the Flow Over an Inclined Flat Plate at Free Stream Mach 
Number 1, NACA TN 3723, August, 1956. 

7 Hadamard, J., Lectures on Cauchy’s Problem in Linear Par- 
tial Differential Equations, Dover Publications, Inc., New York. 

8 Pucci, C., Sui Problemi di Cauchy Non ‘‘Ben Posti,’’ Rendi- 
conti dell’Accademia Nazionale dei Lincei, Classe di Scienze 
Fisiche, Matematiche, Naturali, Serie VIII, Vol. XVIII, No. 5, 
May, 1955. 

9 Zlotnick, M., and Newman, D. J., Theoretical Calculation 
of the Flow on Blunt-Nosed Axisymmetric Bodies in a Hypersonic 
Stream, AVCO Mfg. Co., Rep RAD-TR-2-57-29, September, 
1957. 

10 Garabedian, P. R., Numerical Construction of Detached Shock 
Waves, J. Math. & Phys., Vol. 36, No. 3, pp. 192-205, October, 
1957. 

11 Morawetz, Cathleen S., On the Non-Existence of Continuous 
Transonic Flows Past Profiles, Part I, Comm. Pure and Applied 
Math., Vol. 9, No. 1, pp. 45-68, February, 1956; Part IJ, Comm. 
Pure and Applied Math., Vol. 10, No. 1, pp. 107-131, February, 
1957. 

12,Van Dyke, M. D., The Supersonic Blunt-Body Problem 
Review and Extension, paper presented at the Twenty-Sixth An- 
nual Meeting, IAS, New York, January 27-30, 1958 (available 
as IAS Preprint No. 801). 

13 Sears, W. R., (Ed.), General Theory of High Speed Aerody- 
namics, ‘‘High Speed Aerodynamics and Jet Propulsion,” Vol. VI, 
p. 635, Princeton University Press, Princeton, N. J., 1954. 








The 
ometry 
mount 
shane 
every ] 
a verti 
the gi\ 
area Ol 
less th 

This 
drag } 
corres] 
tributi 
the int 
wing ]| 
optim 
direct] 

The 
Rado1 
The ¢c 
under: 


Cc 
field 
cedur 
The 
gives 
for ] 
conjt 
tain 
ing t 
whic 
1.) 
Wi 
expr 
of it 
conti 


out ¢ 
Amo1 
conti 
Eq. 1 
cal ca 

All 
used 
nish ¢ 
pape! 


and n 





een 
the 
the 
Dar- 
and 
ail, 
ed, 


ga- 


ITO- 


the 
the 
m- 
urd 


ds. 





A Geometric Problem Related to the Optimum 
Distribution of Lift on a Planar Wing in 
Supersonic Flow 


E. W. GRAHAM* 
Douglas Aircraft Company, Inc. 


SUMMARY 


The problem studied may be regarded as a problem of ge- 
Its simplest form (loosely stated) is then as follows: A 


Determine the exact 


ometry. 
mountain rises up from the x-y plane. 
shape of the mountain knowing only the cioss-sectional area of 
every possible cut which can be made through the mountain with 
a vertical plane. In a more complicated version of the problem, 
the given information might be restricted to the cross-sectional 
area of every cut which can be made by a vertical plane inclined 
less than 45° to the y-axis. 

This latter case has direct applications to certain minimum 
drag problems in supersonic flow. The shape of the mountain 
corresponds to the (unknown) shape of the optimum lift dis- 
tribution on a planar wing. The cross-sectional area of a cut is 
the integrated value of the lift along a straight line crossing the 
wing plan form. For a restricted range of line inclinations, these 
optimum integrated lift values can sometimes be determined 
directly. Here it is assumed that they are given. 

The problem in its simplest form was originally solved by 
Radon, who found solutions for a large class of such problems. 
The derivation presented here may perhaps be more readily 


understood 


INTRODUCTION 


ERTAIN minimum drag problems for wings in super- 
C sonic flow can be reduced by the combined-flow- 
field methods of Munk! and Jones*~* (or by other pro- 
cedures) to two-dimensional potential problems.°~!° 
The solution of a two-dimensional potential problem 
gives the minimum drag for a prescribed total lift, or 
for prescribed total volume, etc. It yields also (in 
conjunction with the three-dimensional geometry) cer- 
tain integrated values of lift or thickness, correspond- 
ing to the Hayes’ cuts!! of the optimum distribution, 
which indirectly describe the optimum. (See Appendix 
Ls) 

We will consider planar wings, and derive an explicit 
expression for a lift (or thickness) distribution in terms 
of its known Hayes’ projections and their analytic 


continuation.t We are primarily interested in finding 


Received February 24, 1958. 

* Aerodynamics Research Group 

+ The author is indebted to the JOURNAL reviewers for pointing 
out certain facts which perhaps deserve additional emphasis 
Among these are the following: that the problem of analytic 
continuation has not been completely solved, that H(xo, a 
Eq. 1) cannot be an entirely arbitrary function, and that numeri- 
cal calculations may be quite difficult to make 
Certainly the methods 


(see 


All of these comments seem justified. 
used here are not sufficiently rigorous and comprehensive to fur 
nish answers to all of the questions implied above. However, this 
paper may, as intended, serve asa starting point for more precise 


and more extensive study. 
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optimum distributions of lift from their known projec- 
tions, but the method is not restricted to such cases. 


DEVELOPMENT 


For example, consider the unknown lift distribution 
I(x, y). (See Fig. 1.) The boundary, outside of 
which /(x, vy) = 0, is also unknown. 

Hayes’ cuts are made by planes inclined at the Mach 
angle to the free-stream direction. Our given informa- 
tion is the integrated lift intercepted by any such plane. 
Assuming JJ = 7/2, the traces of these Mach planes 
in the x-y plane are all inclined at 45° or less to the y- 
axis, so a® < 1. Each plane (or its trace) is designated 
by a value of a and x». Our given information is then 


contained in Eq. (1). 

2} . 

| I(x, vy) 6(x — xo — va)dxdy F(x, a (1 
« } wv 


where 6(x — xo — ya) is the Dirac delta function and 


Va. (See Fig. ] for 


so is zero except when x = xo + : 
notation). 

H (xo, a) is known for a? S 
known lift distribution which 
For 1 < a? < ~, H(xo, a) is not known explicitly, but 
is apparently the analytic continuation of the known 
©, FH (xo, a) corresponds to the 


1, and /(x, y) is the un- 


is to be determined. 


function. For a® = 
spanwise lift distribution, which is known from the 
solution of the two-dimensional potential problem. 
The limits of integration, Y and —Y must be large 
enough to cover the full span of the lift distribution. 
Performing the integration in x reduces Eq. (1) to 












} 
f I(xo + ya, y)dy = H(xo, a) (2) 
— VY 
y TAN 
2(xy)=0 
OUTSIDE 
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MACH PLANE 
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Fic. 1. Geometry of the Mach plane cuts 
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It is now convenient to express the unknown func- 
tion /(x, y) as a double Fourier integral as follows: 


(1 ‘Dm? f f fe, w’ eto bw’y) da’ dw (3) 


Replacing x by x» + ya for use in Eq. (2) gives 


l(xo + ya, y) = (1 2m)? f f f(w, w’ ei X 
€ 


pnaara) Ttdey (4) 


I(x, y) 


Inserting this in Eq. (2) and changing the order of 
integration, we get 


ee . - 
(1 an)? | | flo, a'yers f 
= o r —Yy 


dydw'dw = H(xo, a) (5) 


- ee 
1¥(wa-w) 
e- x 


Carrying out the integration in y and inserting the lim- 


its gives 


(1, an)? f f S(, we X 


{2 sin [V(wa + w’)]/(wa + w’)} dw'’dw = H(xo, a) (6) 


The value of Y must be sufficiently large so that the 
entire lift distribution /(x, y) is covered, otherwise VY 
is arbitrary. For convenience, we choose to let Y 
approach ©, It is assumed that we are dealing with 
functions f(w, w’) sufficiently smooth so that 
lim {2 sin [Y(wa + w’)]/(wa + w’)} = 


Yoo 


27 d(wa+w’) (7) 
where 6(wa + w’) behaves as a Dirac delta function and 
so is zero except for w’ = —wa. 

Performing the integration with respect to w’ in Eq. 
(6) by means of Eq. (7), we get!” 


(1/27) f fla, —wa)e™ dw = F(x», a) (S) 
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If a is considered as a parameter, Eq. (8) can be in- 
verted immediately by means of the Fourier Integral 
Theorem to give f(w, —wa) 

f(a, —wa) = | (xo, a)e **** dxo 9) 
Replacing (—wa) by (w’) gives f(w, w’). f(w, w’) can 
then be determined for w’? < w? without the use of ana- 
lytic continuation since /7(x9, a), a? < 1, is given by 
the Hayes’ cuts. 

J(w, +) is also known, from the spanwise load 
distribution, but the unshaded portions of the spec- 
trum (see Fig. 2) must be obtained by analytic con- 
tinuation. 

Assuming that either /7(xo, a) or f(w, w’) can be con- 
tinued analytically, the original integral equation [Eq. 
(1) or Eq. (2)] is essentially solved. Substituting the 
expression for f(w, w’) into Eq. (3), we have 


Mx. 9) = 0 an): f | f 
. e’ WX —w 


ty) dyodw'dw (10) 


Hx», —w’ ‘a) bo 4 


However, it is desirable to reduce the triple integral 
to a double integral, anc, in order to do so, we change 
independent variables and write 


I(x, y) = (1/27)? { f f x 


(xo, a)e?*~ 9% 9%" dxg|w)idwda (11) 


Since /7(xo, a) does not depend on w, the integration in 
w should be performed first. However, the evaluation of 
this integral between infinite limits is awkward because 
of the factor |w!, [see Eq. (12)]. 


is, 9) = a 2m)? f f HT(xo, a) X 
if Gent) olde dxopda (12) 


To avoid this difficulty, we first integrate by parts 
with respect to x», 


f F(x», ae" ~*~ dxy = 


Pe wae Rs ee 


— 1a) xo = 


f oll Oxe” : ‘ re dxo 
: (15) 
= 1) 


For the cases of interest, where /(x, y) is identically 


a 


zero except in a finite region of the x-y plane, 7(+-, 
a) = 0. 
Eq. (12) now becomes 


Ux, y) = [Il in) f [OFT (x0, a)/Ox0] X 
| ” ata ia (leo wd dxpda (14) 


The integral in the brackets can be expressed as 
follows: 





Asst 
smoot rt] 


Sub 


the de 


By 
Dr. A 
ing th 
M. E. 
mum 
which 

Its 
closec 
for a’ 
differs 
plane 
riguez 
of IT 
long | 
lines 
tion ¢ 

If 
nume 


be in- 
tegral 


9) 


) can 
ana- 
n by 


load 
spec- 
con- 


con- 
[Eq. 
the 


10) 


sral 
nge 


in 
| of 


1se 


») 





A GEOMETRIC 


| * i ie odes = (2/2) lim X 


(15) 


‘feos Q (x — xy — ay) — 1]/(x — x — ay)} 


Assuming again that we are dealing with sufficiently 


smooth functions, 


- . 
| ‘ 
| "peal iad | w) de = 


e 
at/ (x — Xo — ay) (16) 


Substituting this result into Eq. (14), we get finally 
the desired expression for /(x, y), 


I(x, v) = (1 27") f f [O17 (x0, a)/OXo| X 


[dxpda/(x — x9 — ay)} (17) 


as the inversion of Eq. (1) (with Y = ©), which we re- 


peat for comparison. 


H(xo, a) = f i) I(x, vy) 6(x% — xe — ay)dxdy (18) 


DISCUSSION 


By direct substitution of Eq. (17) into Eq. (18), 
Dr. A. M. Rodriguez has obtained an identity, verify- 
ing the correctness of Eq. (17). Also, as a check, Dr. 
M. E. Graham has applied Eq. (17) to obtain the opti- 
mum lift distribution on an elliptic plan form wing, 
which is already known from the work of Jones.’ 

It seems probable that whenever //(xo, @) is known in 
closed form for a® < 1, its continuation can be obtained 
for a? > 1. In general, the function /7(a, a) takes on 
different analytic forms in different regions of the xo-a 
plane. Certain examples studied by Graham and Rod- 
riguez suggest that discontinuities in the analytic form 
of I7(xo, a) are propagated unchanged along infinitely 
long straight lines in the xo-a plane. These straight 
lines correspond to singular points in the lift distribu- 
tion on the x-y plane. 

If the known information regarding H(x0, a) 
numerical form, the problem of continuation seems 


is in 
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more difficult since use of the derivatives of J/(xo, a) 
would probably be required. 

After the preceding work was completed, it was 
pointed out by Germain that a large class of such 
problems had previously been solved by Radon! 
whose work is also described by John.'* However, we 
believe that Radon did not consider the problem of 
analytic continuation. The derivation given in this 
paper is, in some respects, less general, but perhaps 
more understandable than the earlier work. 

Alternative methods for using a two-dimensional 
potential solution to obtain the optimum lift distribu- 
tion in approximate form are described and applied in 


references 15, 16, and 17. 


APPENDIX (1) 


The Availability of the Optimum Hayes’ Cuts and Their 
Analytic Continuation 
References 5-10 inclusive give procedures for reduc- 
ing some of the minimum drag problems of three-di- 
mensional supersonic flow to two-dimensional poten- 
The method of reference 7 shows that 
There the 


tial problems. 
the optimum distribution are obtainable. 
optimum combined flow field is constructed using 
simple “‘artificial”’ distributions of singularities. (These 
“artificial” distributions are not confined to the space 
occupied by the wing, and generally do not represent a 
realistic or useful wing in forward flow, hence the term 
“artificial.” ) 

For example, consider a sonic-edge diamond plan 
form and its Mach envelope, as shown in Fig. 3. 

The appropriate artificial distribution of singulari- 
ties shown in Fig. 4 consists of a cylindrical surface of 
sources starting at the envelope rim and extending 
downstream to infinity. (In reference 7 an upstream 
cylindrical source surface of opposite sign is used also, 
but the end result is the same.) The source density is 
constant on any one streamwise element and is given 
by Eq. (1.1) as a function of @ [see Eq. (1-7) of reference 
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Fic. 5. Analytic continuation in the x» — 


7 and note that a sign correction has been made]. 
Equation (1.1) was obtained by solving a two-dimen- 
sional potential problem 


source density = —2} (2/7) cos@ log tan (6/2) — | 

(2/r)+sino!, O56 

source density = —2} (1.1) 
(2/7) cos 6 log itan (6/2) 

(2/r) + sino}, S082 


y 


y 


This singularity distribution produces no disturb- 
ance within the Mach envelope in forward flow. In 
reverse flow a two-dimensional field is created within 
the envelope, since this region is then unaware of the 
termination of the source cylinder at the envelope rim. 
The combined flow field'~‘ is obtained by superimposing 
the perturbation velocities of the forward and reverse 
flow. We add to this combined field the two-dimen- 
sional field created by infinitely long streamwise vortex 
lines corresponding to the spanwise load distribution 
on the wing. The spanwise load distribution is ob- 
tained from the potential jump Ag [see Eq. (I-10) of 
reference 7 |. 


Ag = (2/m) [(1 — y?)/y] log [1 — y)/(1 + y)] (1.2) 


The resultant combined flow field produces constant 
downwash at the surface of the diamond plan form, thus 
satisfying Jones’ conditions? for minimum drag with 
prescribed lift. Also the zone of silence of the diamond 
plan form is undisturbed. It then seems possible that 
this combined field could also have been created by 
some distribution of lift on the diamond plan form. 
It seems sufficient to require that the lift distribution 


should have the same Hayes’ cuts as the artificial dis- 
tribution of sources, which leads to the problem con- 
sidered in this paper. 

The Hayes’ cuts are obtained by cutting with Mach 
planes the source distribution defined in Eq. (1.1). 
(The streamwise vortex lines have no contribution.) 
This gives the function (xo, a) for a? S 1, the shaded 


region of Fig. 5. 


is invariant along that line. This is illustrated in Fig 
5, where the change across line A is H2(xo, a) — Hy(x», 
a). The invariance of this property gives [73(x0, a) — 
[To(xo, a) = Ho(x0, a) — Hy(xo, a) or H3(xo, a) = Hol xo, 
a)— IT (xo, Q). 

At this point a double integration must be carried 
out to obtain /(x, y) from Eq. (17), and it would prob- 
ably be necessary to perform at least one integration 
numerically. This has not yet been done. 

“Artificial’’ distributions of singularities can be used 
i some base-area and volume problems,’ as well as in 
lift problems.’ However, in reference 10, which deals 
directly with the properties of the optimum combined 
flow field, the Hayes’ cuts are derived without recourse 
to any auxiliary singularities. 

A related method has also been used by Jones” to 
obtain a restricted set of Hayes’ cuts which yield the 
chordwise load distribution for the wing. 
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The analytic continuation of this 
particular function into the regions a* > 1 was ob- 
tained by Rodriguez, and is consistent with the tenta- 
tive rule that the change in form of the analytic func- 
tion /7/(xo, a) encountered in crossing any straight line 
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On the Excitation of Pure Natural Modes 
in Aircraft Resonance Testing 


R. W. TRAILL-NASH* 


Australian Defence Scientific Service 


SUMMARY 


\ method is proposed for systematic adjustment of force levels 
to produce a pure undamped natural mode of vibration in multi- 
point excitation of aircraft structures. The question of optimum 
numbers of measuring stations and exciters is examined and some 


conclusions are reached. 


(1) INTRODUCTION 


oun TECHNIQUE of multipoint excitation in aircraft 
resonance testing! offers a means for excitation of 
pure undamped natural modes, but its application 
is hampered by the lack of a systematic method for de- 
termining the appropriate force distributions when the 
damping distribution is not known.j A method is 
proposed here for deriving these force distributions. 
It is quite general and no assumptions as to distribution 
or degree of damping are involved. 

The method requires, initially, the recording of in- 
phase components of displacement at a number of 
points on the structure for each of a number of linearly 
independent force distributions applied at the un- 
damped natural frequency of the mode to be measured. 
From these data a set of linear simultaneous equations 
is set up, the solution of which provides the coefficients 
of a linear combination of the force distributions which 
will serve to excite the pure mode. Arithmetical solu- 
tion of the equations should not be necessary, as a re- 
laxation process can be applied directly to the physical 
system. The method is in essence one of systematic ad- 
justment of force amplitudes, the adjustment being 
based on a set of measured dynamical influence co- 


efficients. 


(2) FORMULATION AND SOLUTION OF THE PROBLEM 


It is assumed that a finite number » degrees of free- 
dom is required for specification of the response of the 
system to an arbitrary distribution of excitation at a 
particular undamped natural frequency w,,.. (This num- 
ber may be considered as the “effective number of de- 
grees of freedom”’ at the frequency w,. It will, in gen- 
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lucidly discussed (see reference 5). 


eral, vary in a particular system with the order m of 
the natural mode.) Displacements, (possibly includ- 
ing rotations) at ” stations on the system are taken as 
generalized coordinates and the (complex) amplitudes 
of these are denoted by the column matrix x. Forces 
(or moments) are applied at » stations and are denoted 
by the column matrix f exp (/w,,¢). It is required to 
establish the force distribution necessary for excitation 
of the pure undamped natural mode corresponding to 
Wm- 

Equations governing the response of the system at 
the undamped natural frequency w,, are written as 


(—wn?A + 1B, + E)x = Th=p (1) 


where A is the generalized inertia matrix; B,, is the 
generalized damping matrix which may contain both 
viscous and structural components (i.e., is of the form 
B,, = w,B + G where B is the viscous component and G 
the structural component), and it is assumed that B,, is 
symmetric; and £ is the generalized stiffness matrix 
and 7'f = pis the set of generalized forces corresponding 
to f. 

For an arbitrary generalized force column the re- 
sponse x is complex and is written as 


yt is (2) 


x= 
and Eq. (1) becomes 
(—w,’A + iB, + E)(y + tz) = p (3) 


In the light of Prof. B. M. Fraeijs de Veubeke’s 
analysis” there is a set of characteristic constant-phase 
linearly independent generalized force columns /, for 
which the corresponding response columns x are con- 
stant phase and linearly independent. One of the 
forced modes is the mth undamped natural mode; 
and this mode is excited in quadrature with the applied 
z is a linear com- 


forces. It follows, therefore that 

bination of the 7 modes, while y is a linear combination 
of (x — 1) modes only, the quadrature mode being ab- 
sent. f 


Suppose that f;, 7 = | . nis a set of linearly inde- 
pendent force distributions which in an experiment are 
applied in turn to the system and x, = y,; + ?2;,7 = 1, 

. n the corresponding measured responses. It is re- 
quired to determine a linear combination of the f; which 


will excite the mth pure undamped natural mode. Such 


t This holds provided that neither all the excitation stations 
nor all the measuring stations lie on a nodal line of any mode. It 
is assumed throughout that this condition obtains 
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a linear combination of force distributions is denoted by 


n 


f= > af; (4) 
j=l 
where the g; are arbitrary multipliers. The correspond- 
ing response is 
n n 


x = Dd ax, = Dd ay; + 22;) (5) 


Jj B 
and Eqs. (4) and (5) are alternatively written, respec- 
tively, as 
f = Fq (6) 
and x = Yq+t 1Zq (7) 


where F is a square matrix of the columns f;, Y is a 
square matrix of the columns y,, Z is a square matrix of 
the columns z;, and g is a column matrix of the multi- 


“j, © 


pliers q;. Eq. (3) then becomes 


(—wn2A + iB, + E)(Vq + iZq) = TFq (8) 


A necessary and sufficient condition for the response 
to be in the pure undamped natural mode is that this 
response be in quadrature with the exciting forces.’ 
This requires that 


Yq = 0 (9) 


Eq. (9) is a set of ” linear homogeneous simultaneous 
equations in x unknowns. The equations have a solu- 
tion, since Y is of rank (m — 1), the columns being linear 
combinations of x — | linearly independent columns as 
noted above. Solution for the g then leads to the re- 


quired force distribution 
f = Fq (10) 


and the corresponding response x = 7Zg, measured 
under this distribution, is the mth undamped natural 
mode. 

Solution of the equations Vg = 0 presents a problem, 
for it is undesirable that such a calculation should be 
required during progress of a test. It may be possible, 
however, to overcome this difficulty by the use of a re- 
laxation technique on the experimental setup itself. 

Suppose that some trial solution g = g, has been 
selected and applied to the system. The real com- 
ponent of response is then 


Yq; = Vs (11) 


and the R/ZS is, in relaxation terminology,’ a set of 
residuals of the equation Yq = 0 corresponding to the 
approximate solution g = g,. The matrix Y consti- 
tutes a set of relaxation operators which indicates the 
increments of g required for the next approximation 
9st+1. These can be applied through the vibrators and 
the new set of residuals y,,; measured. The process can 
be continued until the residuals become acceptably 


small and the required mode is produced. 

Thus, the procedure proposed is initially to deter- 
mine experimentally the columns of Y [Eq. (7)] cor- 
responding to a set of linearly independent force dis- 
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tributions (which in simplest form are single station 
force applications). Systematic adjustment of force 
levels is then carried out experimentally as in conven- 
tional relaxation procedures, the matrix V providing a 
set of relaxation operators upon which successive ad- 
justments are based. When the required condition of 
null in-phase response is achieved, the quadrature com- 
ponent (now total response) is measured. This is the 
desired mode. 


(3) APPLICATION OF THE METHOD 
Choice of Excitation and Measuring Stations 


In practice at the outset the effective number of de- 
grees of freedom of the system is not known, so that 
there is some uncertainty about the number of excita- 
tion and measuring stations required. In order to 
assess the implications of various conditions of excita- 
tion and measurement, the problem is now examined 
in some detail. 

As in Section (2), consider the equations 


(—w,p,"A + 1B, + E)x = p (12) 


and let p;, x;,7 = 1 . n be the set of characteristic 
constant phase linearly independent generalized force 
distributions and corresponding displacement columns, 


Further, put 


x= djk, 7 1k;, 


respectively. 
got... (13) 
where the d; are scalars.* 

Then an arbitrary force distribution p is expressible as 


n 


p= de v;p; = Pv (14) 


where the v; are arbitrary multipliers, v is a column 
matrix of v;, and P is a square matrix of the columns /). 
The corresponding displacement column is 


n n 


x = D>) on; } vj(djk; + ik;) = KDv+ikv (15) 


j=1 j=l 


where A is a square matrix of the columns k; and D is a 
diagonal matrix of the d;. 
Thus, Eq. (12) is alternatively expressed as 


(—wn?A + 1B, + E)(KDv + iKv) = Pv (16 


It is recalled that A and P are nonsingular in view of 
the linear independence of their respective columns, 
while D is a diagonal matrix with one zero diagonal ele- 
ment (that corresponding to the mth mode). 

Consider now excitation at 7 stations. If f denotes 
the applied force distribution, the corresponding gen- 
eralized force column becomes 


bp = Tf (17) 
where 7’ is a matrix of order X r and rank n if n < r, 


and rank rifu>r. Let fs = 1 
linearly independent force columns with corresponding 


r be a set of 


generalized force distributions 


*d; = —cot ¢; of reference 2. 
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EXCITATION OF 


b= 27, £= bb... (18) 


Each of the latter is expressible, as in Eq. (14), as a 
linear combination of the characteristic generalized 


force distributions p;. Thus, 


pb; = P2,, pei...% (19) 


where now the v, are specified columns v. Further, 


? ? 


p= Tf = >> 4.7f: = > ¢.Po, = PVq (20) 
s 1 1 
where the q, are arbitrary multipliers, gis a column ma- 
trix of the g,, isa matrix of the columns 2, and is there- 
fore of order n X rand rank nif nm <rorrankrifn>r. 
This is effectively the introduction of a coordinate 
transformation 


v= Vq (21) 


on the R/ZS of Eq. (16). Eq. (16) can, therefore, be 


rewritten as 
(—w?nA + 1B, + E)(KDVq +1KVq) = PVq_ (22) 


Eq. (22) now allows for an arbitrary number of excita- 
tion stations. It remains to introduce an arbitrary 
number of measuring stations. Suppose that an arbi- 
trary number of stations is selected (possibly including 
those at which the x displacement are measured) and 
that the displacements at these stations are denoted by 


t. It follows then that 
X= Cx (23) 
where C is a matrix of orderh X nand rank nifn <h 


or rank h if n> h. 


expressible as 


Measured displacements are then 


¥ = CKDVq +71CKIq (24) 

It is now possible to examine the outcome of various 

conditions of excitation and measurement, and, for the 

purpose of reference, the relevant sets of quantities are 
restated. These are 


y = CKDIq (25) 
= CKIq (26) 
p = Tf = PV¢ (27) 
and it is required that 
y= 0 (28) 
z=k, ~0 (29) 


The implications of variation in number of measuring 
stations are examined first, and consideration of the 


excitation system follows. 


Consideration of Measuring Stations 


In considering the measuring stations, suppose that 
the excitation system is adequate to provide the desired 
pure mode. It is shown in Section (2) that » exciter 
stations are sufficient for unique solution, and it is 
therefore assumed that this condition obtains. It is 


shown in Section (2) also that either 7 or m — 1 measur- 


PURE NATURAL MODES iii 


ing stations provide sufficient data for solution of the 
problem. The implications of other numbers of meas- 
uring stations are now investigated. 

Number of Measuring Stations Greater Than Number 
of Degrees of Freedom (h>n): ¥ = Oin Eq. (25) isa set 
of homogeneous linear simultaneous equations in which 
the number of equations exceeds the number of co- 
ordinates. However, since C, A, and J’ are of rank n, 


and D is of rank (n — 1), the matrix CAD is of rank 
(x — 1) and the equations are soluble. Further, since 
C is of rank n, omission of any (kh — m) rows |[corres- 


ponding to omission of (# — n) equations in Eq. (25) 
yields 

C, KDVq = 0 
in which C,,, the depleted C matrix, is nonsingular, and 
it follows that 


KDVq =y =0 
which satisfies Eq. (28). In addition, since CA I is of 
rank n, then = ¥ O [Eq. (26)] and s ¥ 0 as required by 
Eq. (29). 
excess of the number of degrees of freedom is satisfac 


Thus, a number of measuring stations in 


tory. 

Number of Measuring Stations Less Than Number of 
Degrees of Freedom (h <n In Eq. (25), ¥ 0 be- 
comes a set of homogeneous linear simultaneous equa- 
tions in which the number of equations is less than the 
number of coordinates. For = (nm — 1) the rank of 
CKDV is (n — 1), and there is a unique solution for the 
g ratios. Solution of CADI’q = O in this case is 
effectively solution of ADIq = 0 with one equation 
omitted, and, since ADV is of rank n — 1, the solution 
exists and is unique [cf. Section (2) ]. 

When hk < (n — 1), CAD V is of rank h, and there is an 
infinity of solutions in which (7 — 1 — h) gq ratios are 
arbitrary. Such solutions do not necessarily satisfy y = 
0. Thus, k < (n — 1) measuring stations are not suf- 
ficient for solution of the problem. 


Consideration of Exciter Stations 


In examination of the exciter question, it is assumed 
that the satisfactory minimum number of measuring 
stations is used; thus, C is taken to be of order 
(n — 1)n and of rank (m — 1). 

It is shown in Section (2) that the required excitation 
can be obtained with the use of as many exciters as 
If the 
number of exciters is r = (7 — 1), CADI is a square 
Therefore, 


there are degrees of freedom of the system. 


matrix of order (7m — 1) and rank (nm — 1) 
y = CKDIq = O has no solution. Similarly, if 7 


“ 


(n — 1) there is no solution. Hence, it is necessary 


that the number of exciters be not less than the number 
of degrees of freedom—.e., r { n. 

Suppose now that r > n. In this case, CAD is of 
rank (x — 1), while the number of coordinates g is 
greater than x. Thus, there is an infinity of solutions 
in which (r — ») coordinates can be arbitrarily speci- 
fied. Such solutions are not necessarily the most ap- 
propriate, as is evident from a consideration of the 
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generalized force column p, Eq. (27). Here there are 


possible conditions 
pb = PVq=0 


fornon-nullg. P is nonsingular so that Il’¢ = 0, which 


makes 


Thus, when the number of exciters is greater than 
the effective number of degrees of freedom, it is possible 
to achieve the condition y = O for a nonzero force dis- 
tribution but accompanied by the condition z = 0—.e., 
a null response resulting from a non-null force distribu- 
tion. In practice, the distribution achieved may consist 
of a component of such a distribution together with a 
component of the required one, resulting in efficient 


use of the forces available. 


Determination of Undamped Natural Frequencies 


The pure undamped natural mode k,, can be excited 
only if excitation is at the corresponding undamped 
natural frequency w,, (and the condition y = 0 can hold 
only for such excitation frequencies). It is important, 
therefore, that the undamped natural frequencies be 
determined with considerable accuracy, and it is sug- 
gested that the method of Kennedy and Pancu‘ is 
suitable. It is conceivable that the first estimate of 
natural frequency, if slightly in error, could be im- 
proved after minimizing the residuals y, of Eq. (11) 
(which has no solution y = O under these conditions) 
and repeating the natural frequency estimation under 
the resulting force distribution. A repetition of the 
process could then yield a more accurate measure of the 
required natural mode. 


(4) CONCLUSIONS 


A technique is proposed for systematic adjustment of 
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force levels to produce a pure undamped natural mode 
of vibration. 

The method requires at least (7 — 1) measuring sta- 
tions when the effective number of degrees of freedom 
is m, but the use of a greater number introduces no 
difficulties. 

Ideally, the number of excitation stations should be 
equal to the effective number of degrees of freedom. 
It is not possible to excite the required mode with fewer, 
and the use of a greater number may lead to inefficient 
use of the forces applied. An extreme case of the latter 
is that in which a completely null response results from 
a non-null force distribution. An initial underestimate 
of the effective number of degrees is therefore desirable 
with subsequent increase until solution becomes possible, 

It is necessary to determine the undamped natural 
frequency with a high degree of accuracy. The method 
of Kennedy and Pancu’‘ is suggested, and there may be 
an advantage in using successive approximations of the 
desired force distribution in improving the frequency 
estimate. 
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On Pitot Pressure in an Almost-Free-Molecule 


Flow—A Physical Theory for Rarefied-Gas 


Flows 


; 


V. C. LIU* 
University of Michigan 


SUMMARY 


\ physical theory of pitot pressure in the transition-flow 
regime—i.e., the moderately-rarefied-gas-flow region—is pro- 
posed. The ratio \/b (the mean free path to the radius of the 
cavity Opening) is assumed to be of the order unity or larger 
4 general formula for the perturbation to the pitot pressure 
in the free-molecule flow is given. This perturbation is at- 
tributed to the intermolecular collisions, which are neglected on 
the basis of the free-molecule hypothesis. The expected rate of 
collision is calculated for rigid spheres, using the classical kinetic 
theory 

Although this is intended as an approximate theory, the 
theoretical results agree surprisingly well with the limited 
experimental data that are available. The present theory shows 
that the ratio Re/\J (Reynolds Number to Mach Number) is 
the governing parameter for determining the intermolecular- 
collision effect on pitot pressure in the transition-flow regime 


SYMBOLS 


c = mean molecular velocity of the effusion beam 
F = approximation to F* 
F* = see Eq. (15) 
WV = Mach Number of the free stream 
\ = rate of molecular flow (number of molecules passing 
per second) 
R = gas constant for unit mass of gas 
Re = Reynolds Number of the free stream 
S(s =m ¢ *"4+4/rs,{1+erf (s;)] 
7 = temperature of the gas 
J = m+C 
J = the relative velocity of two colliding molecules 
V = mean relative velocities of the colliding molecules 
l = radius of orifice 
( = velocity of thermal agitation 
n = molecular number density 
Pp = pressure 
= Bu 


erf(s,) = (2 \ nf, exp (—/*)dt 
0 


= mass velocity of the free-stream molecules along x 
= molecular velocity variable 
= coordinate along axis of pitot tube (measured from 


the plane of the orifice) 


I = see Eq. (A-6) and Fig. 3 
A = (2/3)(1-1/7V 2) 
b = see Eq. (23) and Fig. 3 
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in carrying out the numerical computations 


B = the reciprocal of the most probable velocity 
= ratio of specific heats c,/¢ 
6 = thickness of the wall at the orifice 
€ = |] 745 
bn = mean distance traveled by the free-stream mole 


cules per second 
0 = angle between the velocity vectors of two colliding 


molecules 


K = b/xr 

d = mean free path (following Maxwell 
& = 2x/b 

To = collision cross section 

o = see Eq. (21 

V = see Eq. (10) and Fig. 2 


l = free-stream condition 
F = free-molecule flow 
0 = condition inside the cavity 


(1) INTRODUCTION 


HE ADVENT of high-altitude flight for both ballistic 
G pertenn and sounding rockets for upper-atmosphere 
measurements! has placed a new emphasis on the im- 
portance of the engineering aspects of rarefied gas 
dynamics, which dates back to Maxwell. 

In the case where the mean-free path of the gas, 
\, is much larger than the characteristic dimension, 
b, of the particular flow field in question, the problem 
becomes extremely simplified. The molecules incident 
on, and reflected from, a surface element can be treated 
as having a Maxwellian thermal-velocity distribution 
superimposed upon the mass velocity. This implies 
that the presence of the solid body of size 6 does not 
cause appreciable distortion in the equilibrium-velocity 
distribution of the molecules. This free-molecule 
flow, as it is often called, is readily susceptible to 
analysis, provided the mechanism of molecular re- 
flection at the solid surface is known. Problems of 
free-molecule flow have been thoroughly examined 
(see review by Schaaf*). It has been experimentally 
demonstrated’ that free-molecule flow can be essen- 
tially realized, for special models at least, when the 
value of \/é is still finite. 

Previous work on the development of a theory of 
flow at finite values of \ /b has concentrated mostly on 
the use of special series expansions for the solution to 
Boltzmann's equation for the velocity-distribution 
function.‘~* The calculations are usually formidable, 
and few results of direct aerodynamic interest are 


available. It is our purpose here to propose a physical 
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Fic. 1. Sketch of pitot tube. 


approach to the problem. The basic principle of this 
technique is based on the use of successive approxima- 
tions to bracket the effect of intermolecular collisions. 
This is done without using Boltzmann’s equation.t 
The iteration starts with the free-molecule-flow result 
as its zeroth approximation. From geometrical con- 
siderations one can see that the probability of a re- 
flected molecule colliding with the incident molecules 
two or more times in the neighborhood of the solid 
body is negligibly small compared to the probability 
of its colliding only once; therefore, in the present 
investigation of the almost-free-molecule flows, we 
consider the effect of single collisions between the 
reflected and the incident molecules. 

The analysis of pitot pressure is used as a test case 
for the present theory because (a) the pitot tube has 
been one of the few systems successfully adopted in 
the calibrations of low-density wind tunnels and the 
measurements of upper-atmosphere density from sound- 
ing rockets,’ (b) knowledge concerning the boundary 
conditions of the problem has been reasonably well 
established by molecular-beam experiments,'’ and 
(c) the measurements of pitot pressure in low-density 
tunnels and on sounding rockets may serve as a basis 
for comparison with the theoretical results obtained 
here. 

Precise mathematical theory (such as the yet to be 
obtained solution to Boltzmann’s equation) has its 
place in the final solution of the problem; but it is 
believed that, in dealing with a physical phenomenon, 

+ It should be mentioned here that Heineman’ earlier used 
a similar idea in treating the drag on a flat plate normal to a 
free stream of very high Mach Number. Also, the author’s 
attention has been drawn by one of the reviewers to the work of 
Lune and Lubonski,* who follow the same idea as Heineman’s 


in their analysis of drag problems. 
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a simple conception of the physical nature of the 
process involved is often more conducive to the 
development of new ideas and tests and toward the 
understanding and application of the phenomenon. 

The kinetic concepts concerning free-molecule flow 
and molecular effusion form the root of the present 
theory of an almost-free-molecule flow. This, it is 
hoped, excuses the subsequent rather lengthy and 
detailed discussions of the assumptions involved in 
these simple concepts. 


(2) PiroT PRESSURE IN A FREE-MOLECULE FLOW 


Consider a pitot tube with an orifice of radius 5 
and an axis parallel to the free stream of Mach Number 
M (see Fig. 1). If \1/b > 1, the molecular fluxes in 
and out of the cavity through the orifice are essentially 
independent of each other. 

The rate of molecular flow into the cavity through 
the orifice is given by"! 


MN, = (/4)b?mc fexp (— 512) + YW asi[1 + erf (s1) } (] 


where 7”, and c; denote the molecular density and the 
velocity of thermal agitation of the free-stream mole- 
cules, respectively, and s; represents the ratio of the 
free-stream mass velocity “, and its most probable 
velocity, 1/@1 = (2RT,)'”". 

From the kinetic point of view, the mechanism of 
effusion depends primarily upon the mean free path 
of gas inside the cavity (Ao) relative to the size of the 
orifice. If Xo is of the same order as (or larger than) 
the diameter of the orifice, the spatial and velocity 
distributions of the molecules inside the cavity are 
not sensibly affected by the effusion process. Note 
that \y is usually smaller than A, because of the ram 
effect of the stream. 

The nature of the effusion process is also dependent 
upon the thickness 6 of the walls at the orifice. If 4 
is negligibly small compared to the diameter of the 
orifice, the molecules must effuse independently of 
each other in the form of a molecular stream, each one 
moving with the velocity it had as it came up to the 
orifice. 

Under these conditions,{ the molecular-effusion flux 
from the orifice can be readily calculated to be 


No rb? = NoCo | (2) 


Eq. (2) is found to be in close agreement with measure- 
ments." 
When the steady state of flow is attained at the 


orifice, we have 
Ni = No (3) 
The substitution of Eqs. (1) and (2) into Eq. (3) leads 
to the equation of pitot (or cavity) pressure in a free- 
molecule flow!* 
[(Po pi(7\ 7.) *l, = 
exp (—s17) + Yasill + erf (s:)] (4) 


¢t For molecular effusion from channels of finite length, see 
Clausing’s work.!* 
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PITOT PRESSURE IN AN ALMOST-FREE-MOLECULE FLOW 


Note that 


TIT \V/2 
[moCo/M01|e = [(po/pi)(11/70) "Iz 


and the subscript F denotes ‘‘free-molecule flow.” 


3) GENERAL PRINCIPLE OF THE THEORY OF AN 
ALMOST-FREE-MOLECULE FLOW 


Consider a pitot tube in an almost-free-molecule flow 
(see Fig. 1). The incremental molecular flow generated 
by the collisions between the inlet and outgoing 
molecular fluxes near the orifice (outside the cavity) 
is still small compared to each flux, although not 
negligible. This type of collision will occur, on the 
average, only at a distance from the orifice of the order 
\,;. A consequence of these collisions is the change of 
inlet molecular flux as given by Eq. (1) 

Although an effused molecule may collide again 
with another molecule of the incoming flux before 
it is deflected out of the inlet molecular beam of radius 
b, the probability that this will happen is negligibly 
small compared to the probability of a mere single 
collision. In general, the probability that an effused 
molecule collides Q times near the cavity is obviously 
less than when it collides Q* times, if Q > Q*. The 
present scope of study of the almost-free-molecule 
flows is confined to the single-collision effects. 

In treating problems of an almost-free-molecule 
flow, we take the nominal result of the free-molecule- 
flow theory as its principal value. The deviation from 
this value, contributed by the intermolecular collisions, 
is calculated on the basis of classical kinetic theory. 
The scope of calculation is confined to the collisions 
that occur in a cylindrical region, with the orifice as its 


base, extending upstream to infinity. 


(4) Tue MOTION oF A SINGLE MOLECULE THROUGH 
A GAS 


Maxwellian Distribution 

Consider a single molecule moving at velocity V 
through a gas of density m with Maxwellian distribu- 
tion and zero mass motion. The number of molecules, 
the velocity of which is between v and v + dz, expected 
to collide (per second) with the single molecule (which 
may be called the ‘‘ray’’ molecule, to distinguish it 
from any of the gas with which it collides) is 


4a/ 2o7B,'n, Vv? exp (— 8170") dv (5) 
where ro? = mutual collision cross section; the mean 
relative velocity of the colliding molecules is" 


V, = (1/2) [ V, sin 6d0 (6) 
J 0 


where @ denotes the angle between the velocity vector 
of the “ray’’ molecule and that of the molecule with 
which it collides. We have 

V2? = V2? + v? — QV cos 0 (7) 


After substituting Eq. (7) into Eq. (6), we perform the 


Sl 


integration in Eq. (6) and obtain 
V, = (60V)[(V + 2)? — |V — 2 (S) 
To obtain the total expectation of collisions by a ray 
molecule per unit distance traveled, we integrate and 
then divide Eq. (5) by V. 
E = VV rno’W 


where (see Fig. 2 


exp (—8,°V? ( | ) *Ail 
Vv = : +i2- — exp t*)dt 
B11 B71] J0 


10) 


To apply Eq. (9) in the estimation of collisions made 
by an effused molecule moving through the free-stream 
gas, we let 

V=zamy+C (11 
where “, denotes the free-stream velocity which is in 
the direction of the pitot tube axis, and C, the mean 
velocity of the effused molecules relative to the orifice, 


is given by” 
C= 2 | Bot v4 exp (—Bov*)dv = (3,4)Y 7 Bo (12) 
0 


The distribution function used in Eq. (12) is a modified 
Maxwellian law that holds for molecular beams effused 


from orifices. 


(5) Tue Errect or CAvity PRESSURE ON EFFUSION 


As gas pressure in the cavity increases, an increas- 
ingly large fraction of the effusing molecules collide 
with each other immediately after they emerge. If 
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Collision integrals (function of }b/),). 


Fic. 3. 


we may ignore the fraction of molecules thus collided 
that are deflected into the cavity, it is found that to a 
first approximation the total flux of effusion as given 
by Eq. (2) still holds. The number of effused molecules 
in the cylindrical region for intermolecular collisions 
[see Section (3)] is, however, definitely affected. 

To estimate the effect of the cavity pressure on the 
distribution of effused molecules at various distances 
upstream from the orifice, we proceed with the kinetic 
calculation for rigid spherical molecules which move 
with a uniform velocity cy but in random directions.t 

Let x be the distance from the plane of the orifice 
(see Fig. 1), and N the number of effused molecules 
with velocity co in the beam passing the cross section 
at x per sec. The relative velocity between the 
(effused) beam molecule and a cavity molecule emerging 
at an angle a with the beam molecule is c,? = 2c? X 
(1 — cos a). The number of collisions made by the 
(effused) beam molecules in traveling the distance 
from x to x + dx is!! 


aN = Nao?no(dx a) f C,(1/2) sin ada (13) 
0 
where tan a* = b/x. 


Note that c,dc, = —c? sin a da. 
Eq. (13) 


We obtain from 


N(x)/No = exp [—2o0*mF*(x)] (14) 
where 


F*(x) = (2+/2/3) X 
{1 — (1 + x2/b2)4/[(1 + x2/b2)'/? + x/b]'/7} (15) 


To evaluate NV(x)/No, we use the following approximate 
formula for F*(x): 


F(x) = VY 2A[1 — exp (—e x/b)] (16) 
Vv I 


where A = (2/3)(1 — 1/+/2) and e = 1.745. The 
largest deviation from the correct value is less than 
2 per cent when x/b < 2, within which the most signifi- 
cant contribution by the collisions occur. 

In spite of the crude assumptions used in the analysis 
of scattering at effusion, the results are found to be in 
close agreement with molecular beam measurements 
reported by Istermann, Simpson, and Stern.“ The 
comparison is made for molecular velocity (in beam) 


7 This simplifying assumption is used only in the scattering 


analysis of the effused flux presented here. 
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equal to (3 b)4/ rr, Bo, the value of the mean effusion 
velocity [see Eq. (12)]. 


BETWEEN THE FREE-STREAM AND 
EFFUSED MOLECULES 


(6) COLLISIONS 


Three important factors contribute toward the 


weakening of the effusion-beam intensity. These are 
(1) scattering due to eifused molecules (considered 
above), (2) free (spherical) expansion, and (3) at- 
tenuation by scattering produced by the free-stream 
molecules. Thus, the net beam intensity at &, where 


& = x/b, can be approximately represented by 


[moco/4(1 + &)?] exp [—2o?mbF — (b/A)E] (17 


The loss of total inlet molecular flux in the range 
from x to x + dx due to collisions between the molecules 
of effusion and free stream can be obtained by multiply- 
ing Eq. (17) by Eq. (9). The total loss of inlet mo- 
lecular flow is, therefore, 

n 
Nu = (] 1) rb*nycok [1 (1 +- £)*] yd 
J 0 
exp [—(F/V/2)(b/A1) (m0/m) — (b/A)E|dE (18 


where bn denotes the mean distance traveled by the 
free-stream molecules per second, provided deflection 
of molecules into the cavity from intermolecular 
collisions is ignored [see Section (7)]. To evaluate 
the integral in Eq. (18), we use the approximate 
expression for the exponential function in its integrand, 


namely, 


exp [— (F/V/2)(0/A1) (mo/m)] = 
exp [—A(b/A1)(m0/m1)]} 1 + A(O/M1)(mo/m1) X 
exp [—e(x/b)+ ...]} (19) 


The integration in Eq. (18) can be easily performed 
after the approximation (19) is made. Thus we have, 


after the substitutions \y = (4/2707m) 


and «x = b/\, 


Nu = (1/4) rb nocok exp [—A(b/r1)(M0/m)] X 
[b(x, n) + A(D/d1) (10/1) b(K + €,)}] (20) 


where 


o(<.9) = 1+ (1 + 9) % — x[Ex(—« — ) —- 
Ei(—«)Je“ (21 
For most high-speed flows, 7 > «~!; hence, 


Nu = (1/4) rb nocoE exp [—A(b/A1)(t0/m)] X 
[P(«) + A(D/Ai) (10/11) P(k + €)] (22 


where (see Fig. 3) 


P(x) = 1 + xe" [ [exp (—?) /t]dt (: 


(7) CONTRIBUTION OF INLET MOLECULAR FLOW BY 
INTERMOLECULAR COLLISIONS 


To estimate the net loss of the inlet molecular flow 
with respect to its value based on the free-molecule 
hypothesis as represented by .V; [see Eq. (1)], we must 
subtract from .V;,; the number of molecules per second, 
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emerging from collisions, that are deflected into the 


avity. 
cay xX - 


a 








Assume that the effused molecules move in the mean 
direction x with velocity C. In the case of an adiabatic 
cavity, the quantity C is equal to the mean inlet- & 
molecular velocities. The differential scattering cross 
section in the center-of-mass system,'* which, in the y x 
present case, is at rest relative to the cavity, is o°/4 
neglecting the small-angle forward scattering.'° The 
chance for a scattered molecule being deflected into a ’, -y 
solid angle dQ at an angle 6 relative to its original p 
direction is dQ/4m (see Fig. 4). Note that p . 





dQ = pdpdd cos 6/r* 
Po Spo dd, 


where 
cos8 = —x/[x? + p* + po? — 2ppo cos (¢ — go) eve 
The total number of molecules effused from the Fic. 4. Scattering of the effused molecules. 


orifice that are scattered into the cavity per second is 


. T ‘ . 
N 12 = (*) b nocokt y 4 
| 22x 22 *] 1 n tn* o* xp (— 2n oF — £) 
| | dg | ao | dpx* dp* [ dé ll. acemes = a = — = | (24) 
tr? Jo J0 0 0 J0 (1 + &)? [2 + p**? + po*? — 2p*po* cos (6 — go) |” 


where p* = p/b and po* = po/b. 

Considering the fact that Vy» as given in Eq. (24) will amount to only a few per cent of Ny, in Eq. (22), we may 
facilitate the integration by using a gross approximation to the integrand in Eq. (24). The evaluation of the multiple 
integral in Eq. (24) yields (see Appendix) 


Ny = (2/4)b mocoET (x) exp (—Axno/m1) (25) 


where T is a function of x (see Fig. 3). 
The number of the free-stream molecules emerging from collisions that are deflected into the cavity per second 


is also given by Eq. (25). 
(8) EQUATION FOR THE PITOT PRESSURE IN AN ALMOST-FREE-MOLECULE FLOW 


The condition of equilibrium for mass flow at the orifice is 
No - Nie a= N, —= (Ni = Nie) (26) 


Substituting Eqs. (1), (2), (22), and (25) into Eq. (26), we obtain the equation for pitot pressure in an almost- 
free-molecule flow, for the case 7 >> «~', 
ry ry 1 2 
(po/ pr) (14/70) ] oe 
(Zé) 


[(po/pr)(T1/T 0)" e lL + (xk/V/2m)V (x) exp (—Axno/mi) [P(K) + (M0/m1)AKP(K + €) — 21 (x) 


(Note that [(po/p1)(T1/7»)'’?]» refers to the free-molecule flow.) 

It is noticed that the unknown m now appears on the right-hand side of Eq. (27); 
interpreted only as a governing equation for successive approximations of the pitot pressure. 
tion, we may take the value of mo from the theory of free-molecule flow [see Eq. (4)]. A first approximation to the 


pitot pressure in an almost-free-molecule flow thus becomes, for the case 7 >> x~', 


therefore, Eq. (27) must be 
In the first approxima- 


(po/t1)(T1/T 0)? | 
Po/t . (28) 


[(po/Pr(T1/To) 7 en 1 (x /4/2 9) (x) exp (—AnSA/T1/T 9) [®(«) + (SAY T1/Ty) Ax®(x + €) — 20 (x 


where S = exp (—s517) + V/7s;[1 + erf (s;)]. 


7 SRCRUIERe Aim SeeR ANE the free-stream Mach Number 1/—e.g., at J = 3— 

The ultimate limitation of the present theory is \;/Ao is of the order 5. It is estimated that b/A» = 5 
determined primarily by the size of the mean free path is roughly the upper limit of applicability of Eq. (14). 
of the gas inside the cavity A> relative to the size of the For cases in which / is less than 3, the present ap- 


orifice, 6. In turn, b/Ao is dependent upon b/d; and proximation to the pitot pressure in a rarefied gas 
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Fic. 5. Comparison of pitot pressures in rarefied gases 


should be valid, at least for the range 0 < b/d; < 1. 

It can be easily shown!! that b/\, = 0.499(8/my)'/” 
Re/\J where y denotes the ratio of specific heats of 
the gas; the factor 0.499 is the coefficient in Chapman’s 
formula for viscosity. The functional dependence of 
the pitot pressure upon Re/A/ (or 6/X;) is illustrated 
in Fig. 5. The existing theories of pitot pressure for 
flows at b/A, < 1 and b/d; > 1, respectively, are also 
shown on the same plot. For the sake of comparing 
the present theoretical results with experimental data, 
we use the stagnation temperature of the flow as 7%, 
which denotes the temperature of the gas in the cavity. 
Such an assumption is also adopted in converting the 
results of measurements in low-density wind tun- 
nels.'® ' It is conceivable that with a more-and-more- 
rarefied gas in the cavity, the temperature of the solid 
wall will assume more and more influence in determining 
7). The uncertainty of 75 is the most serious draw- 
back in using existing experimental data (as such) 
for comparison purposes. 

Recent developments in rocket sounding technique 
make it possible to measure the upper-atmospheric 
density below 100 km with fair accuracy.! The 
measurement of pitot pressure as an indirect method of 
obtaining the ambient density at high altitudes has 
been proved satisfactory, as shown by the close agree- 
ment with other established results.!)° Pitot-pressure 
data, measured (near the nose of a sounding rocket’) at 
an altitude around 90 km, would possibly provide a 
valuable check of the present theory for an almost- 
free-molecule flow. Unfortunately, in this sounding 
measurement the gas temperature inside the cavity is 
not available. The apparent variation of the measured 
po/p: with respect to Re/M (or b/d) does exhibit a 
trend which is in qualitative agreement with the 
present theory. 

It is worth noting that at lower speeds of flow, such 
as, say M = 0.512, for which (7)/7»)'/? is approximately 
equal to unity, the theoretical result of the present 
theory seems to show a more consistent agreement 
with measurements,’ '* although at much lower 
speeds, such as J/ = 0.2, Eq. (28) is again ineffective. 
The serious simplification which is introduced by the 
use of Eq. (10) is primarily responsible for this de- 


ficiency. 
Considering the number of simplifying assumptions 
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which were used in the derivation of the final result, 
Eq. (28), which admittedly could be improved by use 
of a more elaborate computation process, we cannot 
expect more than an order-of-magnitude agreement 
with experimental measurements. Yet the agreement 
obtained (see Fig. 5) appears more than encouraging. 
There is, of course, no need to mention that possible 
experimental errors may also tend to aggravate the 
state of uncertainty in checking theoretical with experi- 
mental results. 

Finally, we wish to call attention to the possibly 
general applicability of the present method of approach 
to some problems of rarefied-gas dynamics? thus far 
treated almost exclusively by solving Boltzmann’s 
equation. 


APPENDIX——EVALUATION OF THE COLLISION INTEGRAL 


T*(«) 
- | 2a 2x a1 >] ie 
P(e). = z dd do d po* dp* dé X 
4 7r* 0 J 0 J 0 J 0 J 0 
Ep* po* exp (— 1a0°mF — xé) 
(1 + £) [¢* + pt? + po*? — 2p* py* cos (@ — go) 


(A-1) 


In evaluating the above integral, one should keep 
in mind that T is expected to be of a smaller order of 
magnitude compared to ® at the corresponding value 
of «x. Gross approximations that help to facilitate 
the integration are therefore acceptable. 

Referring to the approximation (18), we take the 


first term of the expansion for exp (—2e°mF). Thus, 


l 2r 2r l l n 
P(e) = , f ago | as | d po* [ dp* [ dé X 
Tr 0 J 0 0 J 0 /7 0 


Ep* po* exp [—Ax(m/m) — xé 
(1 + &)? [E? + p*? + po*? — 2p* po* cos (¢ — do) “ 
(A-2) 


The evaluation of the multiple integral in the se- 
quence ¢, ¢, po*, p*, and & suggests itself where & 
integration may be done graphically if necessary. 

Note that in Eq. (A-2) ¢ and ¢ appear only in (¢ — 
go); we may rotate the coordinate axes so that @o, for 
instance, becomes equal to 7. Then the quantity in 
the integrand involving with ¢ and @» can be written 


[e2? ++ (p* + po*)2]~"”" X 
tp* py* . >|” 
pS - sin (A-3) 
E+ (p* + po*)? = 2 

which is similar to the integrand of elliptical integrals 
of the third kind. If we expand Eq. (A-3) into a power 
series and keep the first term in the expansion, we have 
from Eq. (A-2), after integration with ¢ is performed, 


+ Note added in proof: An application of the present theory 
to the calculation of skin friction on a flat plate has yielded very 


promising results.'9 
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The Effects of Blunt Leading Edges on Delta 
Wings at Mach 5.8+ 


Kenneth F. Nicholson* 

Guggenheim Aeronautical Laboratory, California Institute of 
Technology, Pasadena, Calif. 

July 9, 1958 


— NOTE OUTLINES some of the results obtained by pressure 
measurements on a series of 4 delta wings at M 5.8 in air 
edge Mach Num 


and both sharp and blunt 


Leading edges were subsonic and supersonic 
bers of 0.91 and 1.28 leading-edge 


radius of about 0.5 per cent of root chord Spanwise wing 
sections were triangular to allow some volume for pressure leads, 
with the flat surface of the triangle used as the high-pressure side 
The angle of attack range was from —11.5° to +30° and effects 
of yaw and Reynolds Number were included 

The spanwise pressure distributions turned out to be relatively 
flat and showed no similarity to linear theory. There were no 
large effects of bluntness. At the centerline, the high pressures 
were bounded by Newtonian and exact two-dimensional values 
toward the edge, there were slight changes depending upon 
On the 
low-pressure side of the wing, the pressures tended toward vac- 


Fig. 1 


configuration, angle of attack, and Reynolds Number 


uum. shows a comparison between the experimental 
normal-force coefficients obtained by integrating the pressures 
The 
is attributed to the wedging 
thick 
upper-surface boundary layer, separation occurred only at large 
w/V 


Significant viscous effects appeared in two ways 


and several hypersonic approximations “agreement” 


with linear theory for small w/V 
effect of the triangular section used. 


In spite of a very 


With the 


subsonic leading edge (both sharp and blunt) a chordwise pres 


7 The work discussed in this paper was carried out under the sponsorship 
and with the financial support of the Office of Ordnance Research, U.S. Army 
(Contract No. DA-04-495-Ord-19) 

* Now Aerodynamics Engineer, Douglas Aircraft Co., 


Long Beach, Calif. 
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variation is observed which seems to be caused by an in- 


surt 

crease in local effective span as the boundary layer is carried 
around the edge, combined with pressure gradient and bleed-off 
effects on the boundary-layer thickness. An S-shaped shock-is 
produced by the pressure variations, as shown in Fig. 2. The 
other viscous effect was a uniform reduction in lift with the 
supersonic blunt edges as the Reynolds Number was reduced 


This effect seems to be caused by a “‘viscous pump” action be 


tween the shock and the boundary laver at the edge 


\ll of the viscous-induced pressures were small as compared 
to inviscid pressures producing lift. Some insight into the 
significance of viscous effects at large J can be obtained from the 
two-dimensional form for these effects.?. It turns out that the 


incremental viscous-induced pressure coefficients are independent 
of ./ for large AJ, and depend only on free-stream static condi 
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Fic. 3. Comparison with helium test M 
Here a is free-stream sound speed, L is wing root chord, x is dis 
tance from tne leading edge measured in the free-stream direction, 


and y is free-stream kinematic viscosity 
Pressure data have been obtained previously for delta wings at 
V 2. 


data show a nonconical pressure distribution resulting in pitchup 


13.3 in helium by Bogdonoff and Vas,‘ and their uncorrected 


However, many of their nonconical effects appear to be caused 
by an axial gradient of the tunnel Mach Number \ comparison 
with the interpolated results of reference 3 is shown in Fig. 3 
Differences are attributed to inviscid Mach Number effects rather 
than the difference in the value of the viscous interaction pa 


rameter 1/3/*/ Re 
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Compressible Inviscid Flow Near the End of 


Pointed Afterbodies 


J. W. Reyn 

Department of Aeronautical Engineering, Technological Un 
Delft, Holland 

July 10, 1958 


_. THE SUPERSONIC inviscid flow around a pointed aft- 
erbody is analyzed by means of potential-flow theory, a 


rhe question, 


versity 


singularity is found to exist in the rear point 
whether a shock wave is attached to this point or a stagnation 
Of 


course, the real flow will be determined largely by the boundary 


point is formed, therefore, still remains to be answered 


layer. It is of interest however to consider the inviscid flow in 





788 JOURNAL OF THE AERO 


order to know whether, for instance by means of boundary- 
layer control, the drag of afterbodies could be decreased. For a 
shock wave to originate at the point, a conical-flow solution 
should exist satisfying the boundary condition given by the tan- 
gent cone at the point considered. It can be shown that the 
axisymmetric conical flow around a cone extending upstream to 
infinity does not exist physically because limit lines appear in 
the solution. Axial symmetry is not thought to be essential to 
this result 


pressible-flow region a stagnation point exists at the end of a 


Then it can be concluded that throughout the com 
pointed afterbody. After a brief review of the properties of 
limit lines, a short discussion of the flow with axial symmetry 
around a cone extending upstream will be given. A more de 
tailed discussion may be found in reference 1 


Limit LINES IN COMPRESSIBLE FLows 


A limit line or surface in a compressible flow can be defined as 
a boundary of the flow field, which the streamlines do not cross 
as they are turned back. Such a physically absurd result indi- 
cates that the assamptions underlying the theory are not justi- 
fied. Often in this case, solutions with shock waves can be 
found satisfying the given boundary conditions 

The properties of limit lines have been investigated extensively. 
\ treatise on limit lines in axisymmetric flow was given by 
Tsien.? He finds that for an axisymmetric isentropic flow for 
which a Legrendre transformation is possible this transformation 
becomes singular at the limit line. By this transformation, the 
characteristics in the physical plane are mapped on the charac- 
teristics in the hodograph plane except at the limit line, where the 
characteristics of one family in the physical plane form an enve 
lope, which transforms in the hodograph plane into the limiting 
hodograph. At the limiting hodograph the streamlines are not 
curved and are tangent to the characteristics in the hodograph 
plane. At the limit line the pressure gradient and acceleration 
become infinite whereas the streamline has a cusp. The solution 
cannot be continued beyond the limit line, either with a contin- 


uous or a discontinuous solution 


ConicaAL FLow WiTH AXIAL SYMMETRY AROUND A CONE 
EXTENDING UPSTREAM 


The properties of conical flow with axial symmetry were in- 
vestigated by Busemann* and others. Especially the super- 
sonie flow over the axially symmetric cone at zero angle of attack 
has been given much attention and has become known as the 
Taylor-Maccoll solution. A sketch of the flow field in the 
physical and hodograph plane is given in Fig. 1. It can be 
shown that a limit cone appears in the solution. The flow at the 
limit cone has properties similar to that at the limit line in axisym- 
metric isentropic flows, the results of which, found by Tsien, are 
mentioned previously. A difference however is that Tsien con- 
sidered flow for which a Legrendre transformation is regula 
everywhere, except of course at the limit line, whereas the trans- 
formation of the conical flow is singular throughout the flow field 

The conical flow over a downstream-pointing circular cone ex- 
hibits the same features. From the symmetry of the flow equa- 
tions, it follows that such a flow can be obtained by rotating Fig 











a PHYSICAL PLANE b HODOGRAPH PLANE 


Fic. 1. Conical flow with axial symmetry around a cone extend- 
ing downstream. 
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a PHYSICAL PLANE b HODOGRAPH PLAN 
Fic. 2. Conical flow with axial symmetry around a cone extend- 
ing upstream 
—> 
M> 
Fic. 3. Supersonic inviscid flow around a body of revolution 


pointed at both ends 


l(a) around the y axis and Fig. 1(b) around the uw axis. The result 
is shown in Fig. 2. The flow is now reversed with respect to the 
cone surface. Therefore, and since the velocity normal to the 
radius from the apex of the cone is subsonic in between the cone 
surface and the limit cone, it is impossible to get a shock wave in 
this region. It is therefore concluded that conical flow with axial 
symmetry around a cone extending upstream does not exist. A 
qualitative sketch of the supersonic flow around a body of revo- 
lution pointed at both ends is given in Fig. 3. The Mach Number 
is chosen so large that the bow shock wave is attached to the fore- 
most point. At the tail, a shock wave stands somewhat upstream 
of the rear point, around which a subsonic region exists bounded 


downstream by a sonic line. 
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A Note on One-Dimensional Plasma Motion* 


Robert A. Gross 
Chief Research Engineer, Fairchild Engine Division, 
Fairchild Engine and Airplane Corp., Deer Park, N.Y. 


July 11, 1958 


SYMBOLS 


a speed of sound [see Eq. (8) } 

1 = area 

6 = Alfvén wave speed [see Eq. (8 
B = magnetic-field vector 

c¢ = speed of propagation of small disturbances [see Eq. (8) 
e = internal energy 

E = electric-field vector 

H = magnetic intensity = B/u 
J = current-density vector 
Me= Mach Number = V /< 

p = pressure 


* This work was partially supported by the Air Force Office of Scientific 
Research: under contract No. AF49(638)-15 
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t addition ‘see Eq. 3 
s constant 
m 


ysolute temperature 


ocity in x direction 
V elocity vector 
ensity 


ratio of specific heats 
nagnetic permeability 


ectrical conductivity 


PHYSICS has attracted the increasing attention of 





ntists in many disciplines. Thermonuclear research and 


high-temperature gas dynamics have focused interest on the 


tion of a plasma which is subjected to an external magnetic 
jeld and which exchanges thermal energy with its surroundings 


Because of the complexity of both analysis and experimental 


work in this field it is instructive to examine some relatively 
ple cases to obtain an understanding of some of the effects that 


sim 


expected. This note is written in this vein 


can be 
BASIC Egt ATIONS 


Consider a moving plasma subject to the usual magnetohydro 


dynamic The fluid is considered to be electrically 


assumptions 


neutral (negligible charge separation electromagnetic-wave 


propagation is unimportant 


| negligible displacement currents 
diagonal Maxwell-stress tensor terms are 


significant; 
and the 





the coefficients of viscosity and heat transfer are zero; 


electrical conductivity of the plasma is very large (1/¢ 0 
The equations that govern the motion of such a plasma are 
Conservation of mass: 
Op/dt) + div (pV) = 0 1 
Conservation of momentum: 
6(dV dt) + grad p — (1/u) (curl B) X B = 0 2 


Conservation of energy: 


d/dt) [(pV2/2) + pe + (B2/2y)] 4 
div(pV} (V2/2) + [yRT/(y — 1)]} 4 
1/u) [B2V — (V-B)B]) = p(dQ/dt) (3 


where dQ dt is the rate of heat addition per unit mass 


Maxwell's equations: 


OB/ot) + curl EF = 0 } 

curl H = J 5 

div H = 0) 6 

Ohm’s Law: (J/o) ~E+VXB=0 7 


ONE-DIMENSIONAL STEADY MOTION 


Let us restrict our attention now to the particular case where 


the plasma motion is steady and one dimensional. Thus, 


V = [u(x), 0,0 


In addition, the magnetic field is plane and orthogonal to the 


v direction. Thus, 


We permit thermal energy dQ to enter or leave the plasma by 
means of thermonuclear reaction, radiation, or continued ioni- 
zation or reassociation. For the one-dimensional steady case the 
basic equations become much simpler. The previous equations 


become 
d(pu) = O la 
pudu + dp + (BdB/yp) = 0 2a 
dQ = udu + [yRdT/(y — 1 + (2BdB/up) — (B*dp/yup*) (8a 
Combining Eqs. (4) and (7) gives the frozen-flux relation, 
d(Bu) = 0 da 
Eq. (5) shows that the current flow must be in the y direction 


and of a magnitude given by 


FORUM TSO 


J=-i1 uw dB(x) dx 
Eq. (7) together with the given restrictions on V and B show 
that the electrical field has no sources nor sinks. Thus, 


div E = div (VX B) = 0 


Finally, Eq 6) is satisfied since by definition of our magnetic 


field 
OB.(x)/Oz = 0 ia 


Thus, by confining our attention to conditions defined bv the 
given restrictions on V and B, we automatically satisfy Eqs 


and (6 


It has previously been shown in references 1 and 2 that the 


speed of propagation of small disturbances ini a plasma moving 


in one dimension is given by 
-==aq--+ - Dio + B? up x 


where a is the usual speed of sound and / is the Alfvén wave speed 
Eq. (la), (2a), (8a), and (4a) can now be combined to vield the 


following interesting results 


d(B?) = 2 — 1)BdQ/c?(\f? — 1 a) 
dic? — 1l\)d0 — | x 
2—- ( + l 10 
d(V?) =2 — 1)Af*dQ/(1 — MM 11 
d(M?) = — IAMAdQ/1 — MV » 
; ( 2 - . Wf? — 1); 12 
where i? 2a 


Eqs. (9) through (12) permit us to draw the following conclu 


sions for dQ > O and 1 < 4 < 2. when 
W< 1 heat addition decreases B 


lJ > 1 heat addition increases B 


} > | heat addition increases « 


] l- 2/c7) (2 — J? < 1 heat addition decreases « 
1/5 1 — (b?/c?) (2 — 4 > J? heat addition increases « 
Vf > 1 heat addition decreases u 

VW < 1 heat addition increases u 

WV > 1 heat addition decreases J 

VW 1 heat addition increases \/ 

Eqs. (1) through (12) are identical to the usual aerodynamic 
results if we set 6 = 0. In addition, the conclusions are the 


same as the nonmagnetohydrodynamic case except for the 
effect of heat addition on the speed of propagation of small dis 
turbances. In the subsonic regime, the conditions under which 
heat addition causes ¢ to decrease are related to the magnetic 
field strength If the magnetic field is very strong, such that 
b/c = 1, then in much of the subsonic regime heat addition results 


in a decrease in the speed of small disturbances, « 
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Hypersonic Flow Near the Forward Stagnation 
Point of a Blunt Body of Revolution 


Hakuro Oguchi 

Assistant Professor, Aeronautical Research Institute, 
University of Tokyo, Tokyo, Japan 

July 23, 1958 


N THE PRESENT NOTE we are concerned with the flow near the 


forward stagnation point of a blunt-nosed body of revolution 
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when placed in a hypersonic flow. We are especially concerned 
with an approach to the problem when the thickness of the 
viscous layer is great enough to permit no distinct division of 
the flow field behind the shock wave. 

The fact that, at hypersonic speeds, the ratio k of the density 
ahead of the shock to that behind it is sufficiently small leads 
us to make the following simple assumptions: (1) the flow 
behind the shock wave is incompressible and (2) in the vicinity 
of the stagnation point, r ~ x (r is the distance from the axis 
and x the coordinate along the surface). Let us denote the 
velocity components in the x and y directions (along and normal 
to the surface) by u and v, respectively, then Navier-Stokes 
equations are expressed in terms of vorticity, w[=(Ov/Ox) — 
(Ou/Oy) J, 


u(Ow/Ox) + v(Qw/Oy) (u/x)wo = 


py {(O0%w/Ox2) + [(1/x)(Ow/dx) | — (w/x?) + (0% oy’); (1) 
Setting « = xf’(y),v = —2f(y) so that the continuity equation 


is satisfied and substituting these into Eq. (1), we obtain 
vf” + off’ = 0 (2) 


which becomes by the integration 


of?” + Off” —f'?+ K? =0 (3) 


where A? is a constant of integration 
Let us suppose that the flow just behind the shock wave 


9) fla 


behaves as an inviscid one. Then from Eq. (2 ~ 0 there. 
According to Li and Geigers’ analysis,! the shock-wave conditions 
are given by 

, 2 


f =ku,/2, f’ =u,/Rs ‘ 
grr “ - aty =o ) 
f" = (1 — ku, /kRs? 
where u,, is the undisturbed velocity ahead of the shock wave 
and Rs the radius of curvature of the shock wave, its sign being 
positive when the shock wave is convex to the oncoming flow 
Substituting the above-imposed conditions into Eq. (3), we 
find 
K = Vk(2 — bk) u.,/Rs 
Introducing the transformation 


WV Ky P(n) = f(y), n= V/K vy 

we obtained from Eqs. (2) and (3), respectively, 
6 + 260” = 0 5) 
@’’’ + 266” — o’2?4+ 1 =0 (6) 


The asymptotic solution of Eq. (5) or (6) for large n can be easily 
found to be approximately expressed in the form 


® = a + bn + cn? + O(n) (7) 
where 
1 n’ Par’ ‘zAA sree 
@¢=d exp[—(2an’’’ + bn’’?2 + 
2cn’’’ 3/3) |dn’dn’'dn’"’ 


and a, b, c, and d are the constants of integration. Here it is 





2.0 


, 
pe 

= 

u 

# 








l. = ¢ “a 
. G'— 7G" | 
gp” 1 2s 
0 a | — = 
0 1.0 2.0 2 3.0 


Fic. 1. 


SPACE SCIENCES 


DECEMBER, 1958 





20_ ORs 
18 
$00) | By 


J 5/R;: — 
60) --- 


1.8 























Fic: 3. 


worth noting that the lower bounds of the triple integral may be 
replaced without any significant error by arbitrary but suffi 
ciently large values 

We consider the value of 7, ns, corresponding to the shockwave 
location, sufficiently large so that the integral in Eq. (7) may be 
ns is related 
to das ns = [k(2 — k)]'/44/R.6/Rs where R, is the Reynolds 
Number based on the undisturbed velocity 


ignored, as compared with the quadratic terms. 


ie., Ro. = uRs/v. 
Then, in order to be compatible with the shock conditions (4), 
the solution (7) for large y involving ys must be written as 


®~ AVR, + Bn + Cy2/VR. 8) 


where 
A = [k(2 — k))—/4{(k/2) — (8/Rs) + [11 — &)2/2k] X 
(6/Rs)?} ” 
B = [k(2 — k)]-!/*{1 — [(1 — &)2/k](8/Rs) } 
C = [k(2 — k)]-9/4((1 — &)2/2k] 
Thus we have ®”’ ~ 2C VR, as the condition imposed on 


® at large yn involving ns. For given k and R, the solution of 
Eq. (5) or (6) can be determined so as to satisfy the above condi- 
tion with (0) = &’(0) = 0. 


as the only unknown, is obtained from Eqs. (9) with the value 


The value of ns (or 6/Rs), retained 


of B given, as seen from Eq. (8), by the asymptotic value of 
®’ — n&”’ for large 7. 
enough so that the solution for it may be regarded as the asymp- 


If the value of ys thus obtained is large 


totic one, then the solution may be acceptable as a consistent one 

In the same manner, for the initially given k and ©’’(0), we 
can determine the solution ® and the remaining unknowns R, 
and 6/Rs. (See Fig. 1.) 
for the initial values of k = 0.1, 0.15, and 0.2. 


Actual calculations have been effected 
From the results, 


the values of ’’(0) and 6/Rs are plotted against 1 VR, in 


Fig. 2. Since the local skin friction on the wall, 7,,, is given by 
, g & 

= ’’((0) [k(2 — p)]? 4 1/2 3/3... 

Tw P’’(O) [k(2 ) ]°’ *( pp) (u./Rs)?/ “x 


we can see from this figure that it increases almost linearly with 
1/V R., provided that the other parameters are held fixed. 
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Boundary-Layer Transition With Gas 
injection* 


Cher J. Scott and Gerald E. Anderson 

Rosemount Aeronautical Laboratories, University of Minnesota, 
Minneapolis, Minn 

July 14, 1958 


HE MASS-INJECTION PROCESS has been proposed as a method 
; cooling aerodynamic surfaces, and, since the amount of 
coolant required to maintain practical wall temperatures is con 
siderably larger for turbulent than for laminar boundary layers, 
knowledge of the effect of the cooling method on the transition 
process is certainly important Exploratory studies reported 
here were conducted at Mach Number 3.7 to ascertain the ef 
fects of gas injection on the stability of the laminar boundary 
layer on a conical surface 
The experimental model employed was a porous, 16°, hollow 
cone with a base diameter of 4-in. and a slant height of 14'/s in. 
The conical surface was formed from a sheet of porous, sintered 
stainless steel. According to theory,'» ? the porosity must vary 
as 1/\/x to produce a constant surface temperature with laminar 
flow The cone had a 1-in. solid-steel tip, and, in the first 3 in 
of porous surface, the porosity was low For the remainder of 
the surface, the porosity was within 10 per cent of that specified 
The range of transition positions for the following data is noted 
m the model photograph in Fig. 1. Fig. 1 also shows a photo- 
micrograph of the surface material indicating the distributed na 
ture of the porosity and the roughness. The average roughness 
height as measured by a Brush surface profilometer was 60 


microin., Which is comparable to a production grind finish. The 
plot shown of Stanton Number reduction ratio versus fy, an in 
jection parameter, indicates the cooling possibilities and justifies 
the range of injection rates of the following data. Note that ac 
cording to theory, the Stanton Numbers go to zero at finite 
injection rates. Theory indicates that all gradients normal to 
the wall vanish at this blowing rate, and the boundary layer 
is said to be separated or “blown off”’ from the wall 


2 were obtained from shadow 


The data presented in Fig 
graph photographs and cone-surface temperature measurements. * 
Observed transition Reynolds Numbers, normalized by the tran- 
sition Reynolds Number with no injection and at zero heat 
transfer, which for this case was 3.5 million, are plotted versus 
the equivalent flat-plate injection parameter and also versus 
the area-averaged wall to free-stream temperature ratio ahead of 
transition. A geometrically similar, impermeable, smooth, 
lucite cone produced a transition Reynolds Number of 3.6 


million under the same tunnel operating conditions 
* This research was supported by the USAF under Contract No. AF 18 


600)-1226 monitored by OSR of ARDC. The assistance of D. R. Elgin in 


obtaining these experimental results is gratefully acknowledged 
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The plot on the left in Fig. 2 (data at zero heat transfer) shows 
that while the observed transition Reynolds Numbers do de 
crease with increased injection (more for helium than for 
air), the percentage decreases are not large for injection rates 
which produce significant cooling effects Note that data were 
taken beyond the theoretical laminar blow-off limit. No evi 
dence of blow-off was noted in the shadowgraph photographs for 
these conditions. They did not look characteristically different 
than for the lesser injection cases. The plot on the right re 
sulted from the injection of precooled helium at a constant rate 
of tu — 
caused by injection can be overcome by the cooling possible with 
Note that at a wall to free-stream tem 


0.5 and shows that the Reynolds Number decreases 


this rate of injection 
perature ratio of about two, the transition Reynolds Number is 
greater than at the zero-injection recovery condition. For com- 
parison, the solid-wall infinite stability prediction of van Driest* 
is also plotted 


An extended experimental program is under way 
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Experimental Mixing Profiles of a Mach 2.6 
Free Jet* 


Edward T. Pitkin and Irvin Glassman 

Graduate Assistant and Assistant Professor, Respectively, 
Department of Aeronautical Engineering, Princeton University 
Princeton, N.J 

July 18, 1958 


M* EXPERIMENTAL DATA describing the development of 
subsonic jets and some few data on supersonic diverging 
jets appear in the literature. The only information known to 
the authors on relatively shock-free and parallel-flow supersonic 
jets exhausting into a quiescent gas appears in academic theses.* © 
Experimental pressure and velocity profiles for a Mach 2.6 jet 
are published here since numerous requests for the original theses 
indicate that wider dissemination of these data and results is 
desirable 

The test equipment consisted of a high-pressure storage system 
which provided air at ambient stagnation temperature to an 
axially symmetric Mach 2.6 nozzle of 2.55 in. exit diameter. The 


* This work was sponsored by Project SQUID which is supported by the 
Office of Naval Research under Contract N6-ori-105, T.0 ITI, NR-098-038. 
Reproduction in full or part is permitted for any use of the U.S. Government, 
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nozzle was designed by a characteristic method to give uniform 
parallel flow at the exit plane. Further details of the experi- 
mental apparatus can be found in references 4-6. 

Impact and static pressure traverses were taken with multi- 
probe rakes at increments of ten nozzle-exit radii along the longi- 
tudinal axis of the jet. Impact pressure profiles at the exit plane 
and downstream stations of 0, 20, and 30 nozzle radii are given in 
Figs. 1-3. These profiles suggest the presence of shock waves 
near the jet center line in the region immediately downstream of 
the nozzle. The shocks disappeared further downstream, and no 
effects on the downstream velocity profiles could be attributed to 
them. These observations were verified with more extensive 
data than that shown and with schlieren photographs. Typical 


static pressure profiles also appear in Figs. 1-3. The static-pres 
sure dips in the center were caused by the above mentioned 
shocks 


observed in these profiles outside of the region affected by the 


Considerable variation from ambient pressure can be 
shock waves. The phenomenon is at variance with the usual 
assumption of constant ambient static pressure throughout the 
jet and has been discussed in detail in reference 7. Mean curves 
were drawn through the impact and static-pressure data at each 
station. These values were then used to establish velocity pro- 


files by the Rayleigh pitot tube relationships. The nondimen- 


9 


sional velocity profiles also are given in Figs. 1-3. Essentially 
siinilar velocity profiles are obtained at positions greater than 25 
nozzle radii downstream. 

The boundary of the inner core of full exit veiocity, the locus 
of one-half the centerline velocity, and the centerline velocity 


The core extends 25 radii downstream, 


decay are given in Fig. 4 
considerably longer than the previously reported results for sub- 
sonic jets, and the velocity spreading angle is about 1.5°, consider- 


ably smaller than the previous results.” The initial points on 
this curve show that the shock diamond formation near the exit 
weak and affected the 


Thus, it appears that the mixing of relatively shock-free, parallel 


was relatively velocity only slightly 
flow supersonic jets with quiescent air is a much slower process 
than one would generally anticipate. 

Extensive analyses of these data and others obtained with the 
same equipment and correlations with various theories have been 


presented in reference 5 
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Fic. 4+. Spreading characteristic and centerline velocity decay 
of M = 2.6 jet. R; = nozzle radius, Y = downstream distance, 
R = radial distance, and Ul’; = velocity at nozzle exit 
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The Measurement of Turbulent-Boundary- 
Layer Shear Stress by Means of Surface 
Impact-Pressure Probes* 


Felix W. Fenter and Charles J. Stalmach, Jr. 

Engineering Specialist and Aerodynamics Design Engineer, 
Respectively, Theoretical Aerodynamics, Chance Vought 
Aircraft, Inc., Dallas, Tex 


July 21, 1958 


ee HAS RECENTLY SHOWN that relatively large surface 


pitot probes which extend well outside the viscous 
sublayer can be used to indicate turbulent shear stress in in- 
compressible pipe flow. Preston made use of the fact that the 
law of the wall appears to be universally valid in the neighborhood 


of a smooth surface. The law of the wall may be expressed as 


UV p r=fl(y v V 7/p] (1) 
The function f is assumed to be accurately described by the 
numerical values tabulated by Coles.2 As an extension to Pres 
ton’s work with pipeflow, Hsu* has demonstrated, both experi 
mentally and theoretically, that surface-impact probes can also 
be used to measure the local turbulent skin friction of external 
boundary layers. Furthermore, Hsu successfully confirmed 
Preston's observation that the effect of pitot-tube wall thickness 
was negligibly small 
The surface-impact-probe technique was recently extended to 
include the effects of compressibility’ by assuming the validity of 


* This work was performed while the authors were associated with the 
Defense Research Laboratory of the University of Texas. The research 
was supported by the U.S. Navy Bureau of Ordnance under Contract 
N Ord -16498 
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a law of the wall for compressible Yow over insulated surfaces, as 


follows: 


(U; V/0)V op, Vv ¢ U/l = tly Veo Vr Du (2 


ao = {(4 1)/2) Mi=/41 + [ 1)/2] 3N?} 


7 sin 
where 


and where the subscripts 1 and 7 indicate free-stream and wall 


conditions, respectively. The form of the above law-of-the-wall 


equation results from a mixing-length analysis which assumes 


isoenergetic flow. The function / in Eq. (2) is the same function 
which appears in Eq. (1 

The present investigation was limited to round impact probes 
which lie against the surface The pressure indicated by such a 
probe can be approximated by integrating the local impact pres 
sures over the area of the probe opening. Since the results of 
Preston and Hsu indicate that the effect of wall thickness is 
negligible, the limiting case of b/d = O is considered first, where 
b and d are the inside and outside diameters of the probe, re 
spectively. For this case, the velocity, U’,, indicated by the probe 
reading, corresponds to the boundary-layer mean velocity at a 
distance d/2 from the wall. The relationship between this probe 
measurement and the local shear stress can now be obtained 
directly from Eq. (2) by assuming that the flow is isoenergetic 
The resulting relationship between the surface impact probe 
measurement and the local skin-friction coefficient may be 


written, as follows, for the case of b/d = 0: 


f(Mn)/VWo] R, sin V o( U,/l = V fox Mi)Rj2C,;/2 X 
t[V ful M,)R2C,/8} (3 
where f(M,) = $1 + [(, 1)/2] rM,2} wt1) 
f(M,) = 314 5 1)/2|rM,2} — Cet) 


and R, is the Reynolds Number based on probe diameter. In 
these equations, 7 is the temperature recovery factor and w 
is the exponent in the viscosity power law, u,,/m. = (7,,/7T))* 

The other limiting case of b/d = 1 was also investigated by 
integrating the local imipact pressures over the opening in the 


surface probe. It was found that this laborious procedure 


9 


yielded results which were negligibly different from Eq. (3), as 


long as 


f 
fi(.M,)( Ra/V oa) < 104 
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Fic. 2. A summary of surface-impact-probe data and direct 
skin-friction measurements for various probe diameters, Mach 
Numbers, and Reynolds Numbers. 


This result is consistent with the findings of Preston and Hsu. 

The accuracy of the law of the wall for compressible flow, as 
expressed by Eq. (2), was investigated experimentally in the 
Defense Research Laboratory Variable Mach Number Wind 
Tunnel. The investigation included impact-probe surveys of 
the bouudary layer on the tunnel wall and direct measurements of 
the local shear stress by means of skin-friction balances. Typical 
results are shown in Fig. 1, and it can be seen that good agreement 
between Eq. (2) and the data was obtained in the neighborhood 
of the wall for all values of Mach Number and Reynolds Number, 
based on the boundary-layer momentum thickness. Surface- 
impact-probe meastrements and direct shear-stress measure- 
ments were made at nine Mach Numbers and with four probe 
diameters. Several Reynolds Numbers were utilized for each 
combination of Mach Number me probe diameter. The result- 
ing data are compared with Eq. (3) in Fig. 2 

It is to be concluded that ae surface-impact-probe technique 
can be used to measure with good accuracy the local skin friction 
in supersonic flow over thermally insulated surfaces. For the 
case of isobaric surfaces, the maximum probe size that can be 
used in a given application has been theoretically calculated and 
experimentally confirmed in reference 4. 
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Laminar Heat-Transfer and Pressure 
Measurements Over Blunt-Nosed Cones at 
Large Angle of Attack* 


Victor Zakkay 

Research Associate, Polytechnic Institute of Brooklyn 
Freeport, N.Y 

July 24, 1958 


Bom HAVE BEEN CONDUCTED at a Mach Number of 6, in the 
PIBAL Hypersonic Facility, in order to determine the heat 


* This research was conducted under the sponsorship of contract AF 
33(616)-3265 administered by the Aeronautical Research Laboratory, 
Wright Air Development Center, USAF, under the guidance of Dr. Antonio 
Ferri. 
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transfer and pressure distributions over a slender blunted cone 
at angles of attack of 5, 10, and 15° 
in the tunnel were 615 psia and approximately 1700°R. for all 
tests. The temperature ratio, stagnation to wall, was approxi 
mately 2.8. The Reynolds Number per inch based on fre¢ 


The stagnation conditions 


stream conditions in the test was 0.31 &K 108 The model tested 
has a spherical nose diameter of 1.0 in., a base diameter of 3.75 
in., and a cone half-angle of 20°. The measurements wer 
made at 5 peripheral stations on the model at values of Yy = 0°, 
45°, 90°, 185°, and 180°; y = 180° locates the windward sid« 

In this note the experimental results at a 15° angle of attack 
are presented. A more detailed analysis of the results for all 
angles of attack is presented in reference 1. 

In Fig. 1, the surface static-pressure distributions for the 
various peripheral stations are presented. The pressure dis- 
tributions over the spherical part of the body agree favorably 
with the Newtonian cos? # variation up to # = 55°, and with 
the Prandtl-Meyer expansion theory from 3 = 55° to the tang 
ency point (% being the radial angle measured from the stagna 
tion point). At y = 180° the sphere-cone tangency point is 
55° away from the stagnation point and the Newtonian theory 
yields P/Po = 0.329. For y = 0 this position is 85° away from 
the stagnation point and the Newtonian plus Prandtl-Meyer 
expansion theory yields P/Po = 0.0623. Both values agree 
favorably with the measured values, as observed from Fig. 1. 
On the conical skirt, the flow overexpands with respect to the 
Kopal cone pressure and then is recompressed to approximately 
the Kopal value. At yy = 90° the pressure distribution is ob- 
served to have the same values as the zero angle of attack case 
In Fig. 2 the peripheral pressure distribution at S/r = 7.9 is 
compared with the values obtained by Kopal for cones at large 
angle of attack. The peripheral pressure distributions for a 
given angle of attack can be fairly represented by the second-order 
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ipproximation given by the following equation: 
P=P+Acosy B+ Ceos2y I 


A, B, and 


C are constants evaluated from the pressure measurement at yy = 


P being the pressure at zero angle-of-attack condition; 


2 


(90°, and 180 Included in Fig. 2 is the peripheral pressure 


distribution obtained with the above method. 
rhe heat-transfer results are presented in Fig. 3 in the form of 


the square root of the Reynolds 


Nusselt Number 
Number defined as 


livided by 


where 


model nose radius 
local heat flux 
h,, = stagnation enthalpy 
= stagnation conditions 


= wallenthalpy 


\n estimate of the laminar heat-transfer distribution over the 
body at the zero angle-of-attack position can be obtained by tak 


reference 2, g being the Jocal heat 


ing the distribution of 
flux in B.t.u./ft.2 sec., and gq, is the stagnation point heat transfer 
obtained from reference 3. The heat-transfer distribution ob 


tained in this manner agrees favorably with the experimental 
points. The 


is 10 per cent lower than the experimental value shown in Fig. 3 


theoretical value of g, obtained from reference 3 


\t wy = 180°, an estimate of the heat-transfer rates assuming 
zero cross flow was made using the pressure distribution obtained 
it-transfer values 40 per 


The theory 


from experiments. This resulted in h¢ 


cent lower than the experimental values of refer 
ence 4 for cones at very large angle of attack was then applied at 
y 


heat-transfer results 18 per 


180°, using the local measured pressure values. This yielded 


cent lower than the experiment 
values 

It appears from the above results that for engineering purposes 
a satisfactory rapid estimate of the heat transfer can be made 
at Y = 180° by the method of reference 1, using the local meas 


ured pressure values. 
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Inelastic Stability Theory for Creep Buckling 
of Plates and Shells Under Transient 
Loading 


John M. Bozajian 

Technical Staff, Systems Development Laboratorie 
Hughes Aircraft Co., Culver City, Calif 

August 4, 1958 


and shell structures i 
local 


achieved when ge 


many thin plate 


"ye EP BUCKLING in 
characterized by the 


This behavior 


occurrence of and relatively 


sudden buckling is usually 


metrical and material imperfections are small, of short wave 


length, and randomly distributed. Gerard! has proposed an 


approximate classical stability approach in which solutions for 


1] 


plastic buckling of plates and shells?~* are extended to creep 


buckling through the’ introduction of a time-dependent strain 
The 


placement of the short-time tangent modulus with the dynami 


intensity purpose of this note is to suggest: a) a re 


elastic modulus which is consistent with the assumptions of 


sudden incremental deformation at buckling and (particular] 


for higher temperatures) time-dependency of inelastic deforma 


tion and (b) a generalization in the definition of the secant 
modulus which permits estimates of shell creep buckling under 
transient loads and temperatures. The present remarks have 
evolved as an extension to general shell buckling of correlations 
of recent tests of cylindrical shells conducted at the Hughes 


\ircraft Co.' 


In the stability concept of Gerard! and Bozajian and Tsao’ 

the creep buckling process is idealized as a relatively long period 
deformation of the midplane of the shell (negligible 
followed by 


Gerard applied a time-dependent secant modulus, £ 
stability 


of creep 


bending) i very rapid* buckling period at the critical 


time , and 


short-time tangent modulus, /;, to the existing plastic 


solutions. The tangent modulus was used to describe the incre 


mental behavior during buckling. Gerard notes, with constant 


load during the buckling period, the existence also of some elastic 


unloading which would indicate the use of a reduced modulus, 


but he retains the tangent modulus for the purpose of unifying 
creep and plastic buckling 


It is the writer’s hypothesis that unification might be achieved 


for thin plate and shell buckling in an alternate manner rhe 


major (if not all) inelastic deformations in most materials are 


time dependent, particularly at higher temperatures, and incre 


mental deformations, if sufficiently rapid,* approach elastic 
behavior (linear, instantaneous, and reversible) at the tmception 
of buckling. Creep data for QQ-M-44H magnesium‘ show 


time dependence of the inelastic strains even at room temper 


ature for times measured from a few seconds to several minutes 


For many materials, usual short-time stress-strain tests con 


ducted at constant strain (or stress) rate in a test machine are not 


‘short-time’ enough and thus contain the accrual of inelastic 


creep) strains \t higher temperatures, the creep deformation 


under a low strain or stress rate loading may also obscure the d¢ 


termination of the true elastic modulus These remarks may 


be coupled to a stability concept for inelastic shell buckling which 
} I & 


implies (in the limit) discontinuous behavior at buckling (é€ — 


+ Unification is achieved by regarding short-time plastic 


buckling as transient creep buckling (under variable stress) and 


the consideration of inelastic buckling as, basically, a_ time 


dependent phenomenon 

The utilization of the above model for inelastic buckling leads 
to simple formal procedures for the adaptation of previous plastic 
stability solutions derived by Stowell,? Bijlaard,* Gerard, nd 


others \s proposed by Gerard! the secant modulus is 


stress and strain in 


tensities defined in the deformation theory of plasticity Th 


where for combined stresses o and ¢ are the 


€ 


* This means rapid in comparison to creep rate but slow enough to pre 


clude inertial effects 
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strain is composed of the sum of an elastic, o/E, and inelastic 
(time-dependent), ¢;(t), part. E is Young’s (dynamic) modulus 
of elasticity which is in general a function of temperature. The 
strain rate for many materials may be expressed, negle.. ng 
thermal strains, as 


é = [d(o/E) /dt| + & (a, «€,) (2) 


Analysis?* of the incremental variations in stress and strain 
during the buckling process yields, in the present notation, the 
quantities o/e and do/de. The forimer is defined as the secant 
modulus in plasticity and acquires time dependence, Eq. (1) for 
creep buckling. The latter is the tangent modulus (for de > 0) 
in plasticity and for creep buckling, using Eq. (2), is 


do/de = E(1 — { [ei(o, «:) — (0/E) E/E}/e) (3) 
Since, by the stability concept, the buckling increments are 
and, therefore, 


lim do/de=E 


€— + 


sudden, €— + 


The formal rule evolves as a substitution of the dynamic elastic 
modulus, £, for the tangent modulus, /#;, along with the incor 
poration of the time-dependent secant modulus, £,(¢), in the plxstic 
stability results 

\s an example, Gerard” ® has derived an equation for buckling 


of a long cylinder in axial compression, 


o./(E,.6/R) = 0.6[(1 y.*)/(1 - pt) J? (Ee E,,.)}/? (4) 


=() 


where 6 and R are the shell thickness and radius, v, (==0.3) is the 


elastic Poisson’s ratio, and v the engineering value for the inelastic 


range (0.8 < vy < 0.5). The subscript ¢ denotes the value at 


buckling. For axial compression, the stress intensity (defining 
E,) is the applied axial stress o 3y the approximate theory 
presented here, E,is replaced by E 

Bozajian and Tsao® 7 have reported tests of QO-M-44H mag 
nesium alloy circular cylinders in axial compression at constant 
load and with load increasing linearly at various rates. Temper 
atures ranged from 75°F. to 60°F. and buckling times from 5 
sec. to 5 min. These data were correlated statistically in the 


form 
o./E,(B y 6/R) = a(E/E,,)” a) 
where BE, = oft) /et 6) 


is a generalized secant modulus applving when stress varies with 
time. 8 is a failure probability factor and equals unity for the 
statistical mean value. y is a factor (y = 1 for R/6 = 88) ob 
tained from elastic buckling data that allows for some additional 
dependence on 6/R (probably due to scaling the effects of imper 
fections). In reducing data for variable stress, e(f) wes deter 


mined analytically from Eq. (7) derived as a generalization for 


variable stress of the empirical, constant stress, transient (pri 
mary) creep law given by letting o = constant. 
t K 
e(t) = (o/E) + f(a) /* dt 7) 
J 
Good agreement was found in the appropriate stress and time 
ranges by letting f(a) = bo", (n > 1), and0 < AK <1 with » 
decreasing and A increasing with temperature. The valve of 
E used was a dynamic modulus measured by vibrating canti 
levered strips. 
The values for a@ and p were 0.279 and 0.513 in reference 6 
which reported tests of 56 shells at R/6 = 8). 
additional shells at R/6 = 40 were included in reference 7 and the 


Results for eleven 


combined correlation of 67 shells yielded a = 0.278 and p = 
0.502. Reference 7 also contains a correlation of 59 shells tested 
in pure bending. The functional forms for buckling in axial 
compression given by Eq. (4) (with E,,-—~ E) and Eq. (5) are in 
close agreement thus lending experimental support to the theo 
retical approach. The use of the generalized secant modulus is 
therefore recommended in the estimation of creep buckling 


under variable loads for other shells. The lack of agreement 


between the coefficients, roughly 0.6 in Eq. (4) and 0.3 in Eq. 
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(5), reflects the well-known inability of classical elastic stabilit 
theory to yield the proper value due to the effects of large ¢ 
flections and initial imperfections. It is suggested for shells, jy 
general, that suitable corrections for the elastic case be incorp 
rated if indicated either from advanced theory or experiment 
when applying the inelastic modifications. Experimental dat 
are also invaluable in providing the statistical variations r 
quired in the pursuit of a quantitative index of structural relj 
ability to replace the margin of ‘‘safety Graphical solutions 
have been developed for creep buckling of cylindrical shells under 
variable stress and temperature using isochronous stress-strain 
curves and the concept of stability boundaries.£:7 The sam 
methods may easily be applied to other shell buckling problems 
In closing, it is noted that the stability approach should apply 
only to thin plates and shells exhibiting quasi-stable modes of 
failure in the sense that overall structural configuration is pre 
served, imperfections are small and random, the principal lateral 
deformations are local and sudden as distinguished from a more 
gradual, large-scale primary collapse The stability concept 
should not apply, therefore, to primary creep buckling of long, 
unsupported plates or columns. Indeed, the present results 
would indicate buckling stress proportional to E while Gerard’s 
method! would show dependence on the short-time tangent 
modulus. Neither result provides a proper influence of time on 
inelastic column buckling within the framework of an approximat 


Stability’ theory. 
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A Modified Model of the Supersonic Laminar 
Wake 
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gaan IGATIONS OF THE SUPERSONIC base pressure in blunt-body 
wakes at low Reynolds Numbers have pointed out an appar- 
ent discrepancy between the axisymmetric and the two-dimen 
sional case. According to Chapman's! approximate theory, the 
base pressure in totally laminar two-dimensional flow (when 
transition occurs downstream from the “critical” point of closure 
of the wake) should remain independent of Reynolds Number 
Crocco and Lees,‘ on the contrary, predict an “increased lami 
nar mixing rate’’ with decreasing Reynolds Number, which leads 


to a decrease in base pressure. Experiments in two-dimensional 


flow appear to verify the first result!) ? while the investigation 
of cone cylinders’ indicates a decrease of the base pressure in 
accordance with the second hypothesis. Both theories are based 

* The author would like to acknowledge the contribution of R. Roberson 
to the experiments reported here 
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observations available to date 
The 


laver, neglected in both 


on t dimensional flow, while 
seem to associate the phenomenon with axisymmetric flow 
thickness of the oncoming boundary 
analyses, does not explain the trend It has been found that. it 
increases the base pressure 

It would seem necessary to modify the approximate models of 
the recompression process and to postulate a mechanism which 
can explain the transition from the constant laminar base-pres 
sure ‘“‘plateau’’ to the decreasing base pressure at low Reynolds 
Numbers 
internal flow in the closed 


As viscous effects become dominant, the dissipation 


in the wake increases. The mass 


returned from the recompression point, which is the mass ‘“‘scay 


enged”’ by the shear laver, increases as Re The return-flow 
velocity must increase accordingly The dissipation per unit 
breadth of a two-dimensional wake in the channel-like return 


flow varies qualitatively as (see sketch on Fig. 1 for a definition 


of the svmbols 


o u(du;/dy®Lh ~ pa /hPLh ~ plr(1/R Lx?/h l 


where the mass scavenged by the shear laver is (7;/), and it is as 


sumed that the similar asymptotic laminar shear-laver profile 


still holds. This can be viewed as a minimum which can in 
crease by internal turbulence, if necessary, to dispose of the 
acquired by the wake, but it cannot decrease rhe work 


on the wake by the traction at the separating streamline 





in the shear laver is (the two-dimensional wake is assumed to be 


isobaric, which is well substantiated by experiments 


d(u/l 
Tux ~ pu? \ R 1 + 2 
d V/XN R é 0 


reduced velocity gradient and the reduced velocity 


ire known from Chapman's solution. Equating, it follows 


with Ll ~a 


pproximately, 


ion angle 





The length of the shear laver decreases and the e: 


increases with decreasing Revnolds Number, resulting in a lower 
But then tl 


base pressure, as observed recompression dynam 





ics postulated by Chapman cannot hold any longer; that i 
single trailing shock causes more pressure rise than 1s required 
to prevent the escape of mass scavenged by the shear laver from 


] 


the wak« However, if the recompression occurred through the 


intermediary of several shocks and the wake consisted of several 
cells as depicted on the sketch, each just stopping the separating 


streamline, then a direct modification of Chapman's analysis 


gives for the overall base pressure 


where » is the number of discrete recompression steps, and the 


development of the shear layer in any 


interval still follows the 
asymptotic solution 

Such a ‘cellular’? wake would satisfy simultaneously the con 
dition of conservation of energy Eq. (3) and the conservation of 
mass in the wake. As the Reynolds Number decreases to the 
point where internal dissipation prevents the return of the sca 


venged mass in a single-shock wake, an increasing number of cells 








appear with the attendant decrease in base pressure as given by 
Eq. (4 
recompression region widens, and the shear-layer profiles are 


In reality, the process is obviously continuous, the 
accordingly affected; however, the multiple-cell recompression 
model appears to supply an adequate phenomenological picture 

Fig. 1 shows an extrapolation of Kavanau’s experiments by 
the qualitative rule Eq. (3) with the single point at Re = 10‘ used 
to evaluate the proportionality constant Assuming a single- 
shock structure at the ‘‘plateau”’ (and evaluating the axisymme- 
(4) at that point), the required 
This 


suggests that the decrease in base pressure is steeper than indi- 


tric mixing constant ux by Eq 


number of hypothetical ‘‘cells’” is shown on the curve 


cated by Kavanau, and that it begins at a lower Reynolds 


Number 
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Some measurements on 4:1 cones are also show Unfortu 
nately, very small models had to be used to obtain these Revnolds 
Numbers rhe pressure was measured through the ting, 

in. length of 0.015-in. tubing followed b l-in. length of 
O.032-1n tubing ind then a conical hold Per ps the gradual 


increase in diameter and the fact that the sting length to base 





diameter on models 2 and 38 differed without differences in the 
measured base pressure can be taken as an indication that sting 
interference effects were absent Separate leakage and res] c 
time tests on the assembly showed that the errors caused by the 
pressure-transfer system were within 50 mict below the least 
count of the measurement 
rhe experimental points fall below Kavanau’s data because of 
difference in model geometr\ Chapmaz malvsis of the 
pote tial isvinmetric flow see Table 2 and 3 of reference 5 
poiits out that the “equivalent” free-stream Mach Number for 
L con 0° half-angle) is 4 per cent higher and the equi\ lent 
pressure IS per cent lower than free strean Combining these 
correctious Wit the basic base pressure relation ur data can 
be brought to agree with the cone evlinder 
The trend indicated by these experiments is consistent in show 
radual slope ilmost a constant plateau. Although 
icv cannot be claimed for these tests, one may also 
infer that the decrease in base pressure occurs at lower Reynolds 
Numbers than previously inferred. However, the fact that the 


base pressure decreases, as shown by Kavanau, seems real and 


justifiable It is not a characteristic to axisymmetric 
flows, alth 


it a higher Revnolds Number in this case than in two dimensions 


peculiar 
ugh one would expect the end of the plateau to occur 
This is because of the area contraction of a conical wake toward 


the recompression regioa, higher internal gradients, and larger 


dissipation. Also the gradual recompression of the conical flow, 


as compared with the discontinuous plain recompression shock, is 


likely to permit more gradual transition from the plateau to 


the dissipation-controlled regime 
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